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ABSTRACT 

We present the first study to synthesize results from five different exoplanet surveys using three inde¬ 
pendent detection methods: microlensing, radial velocity, and direct imaging. The constraints derived 
herein represent the most comprehensive picture of the demographics of large-separation (> 2 AU) 
planets orbiting the most common stars in our Galaxy that has been constructed to date. We as¬ 
sume a simple, joint power-law planet distribution function of the form d 2 N v \/(d\ogm p dloga) = 
A(m p /Msat)“(a/2.5 AU)^ with an outer cutoff radius of the separation distribution function of a ou t- 
Generating populations of planets from these models and mapping them into the relevant observables 
for each survey, we use actual or estimated detection sensitivities to determine the expected observa¬ 
tions for each survey. Comparing with the reported results, we derive constraints on the parameters 
{a, /3, A, aw} that describe a single population of planets that is simultaneously consistent with the 
results of microlensing, RV, and direct imaging surveys. We find median and 68% confindence inter¬ 
vals of a = -0.86to;i9 (-0.85±g;ig), /3 = l.lj^ (l.l±J-j), A = 0.21±g;?§ dex" 2 (0.21±g;f§ dex" 2 ), 
and a ou t = IO-47 AU (12+g° 2 AU) assuming “hot-start” (“cold-start”) planet evolutionary models. 

These values are consistent with all current knowledge of planets on orbits beyond ~ 2 AU around M 
dwarfs. 

Keywords: methods: statistical - planets and satellites: general - gravitational lensing: micro - 
techniques: radial velocities - techniques: high angular resolution - stars: low-mass 


1. INTRODUCTION 

Understanding the demographics of exoplanets is req¬ 
uisite to the development of observationally-constrained 
formation and migration models. Exoplanet discovery 
surveys have revealed a large diversity of systems, many 
of which look nothing like our own Solar system. The 
overwhelming majority of planets discovered through 
such surveys have been detected indirectly by monitoring 
variations in: 1) the centroids of absorption lines in host 
star spectra (RV), 2) the apparent brightness of the host 
star (transit), or 3) the magnification of a background 
source near our line of sight to the host star (microlens¬ 
ing). The physical processes that produce these observ¬ 
ables are fundamentally different, and consequently, the 
properties of the planetary systems they reveal are also 
different. In other words, each technique is sensitive to 
planets in a given region of parameter space (e.g. mass- 
orbital period space), with varying degrees of overlap 
with regions probed by other techniques. 

In general, RV surveys are most sensitive to short- 
period, massive planets, transit surveys are most sensi¬ 
tive to very short-period planets with large radii, and 
microlensing surveys are most sensitive to a range of 
planet masses and orbital separations near the ice line 
(~ few AU). Practical limitations currently restrict the 
sensitivity of direct imaging searches to the most massive 
(> Mj up ) planets at large separations (> 10 AU). Coro- 
nagraphy and high angular resolution obtained through 
the use of adaptive optics are required to reveal plan¬ 
etary companions that are several orders-of-magnitude 
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fainter than, and located at small angular separations 
from, their hosts. Although planet detections from in¬ 
dividual methods are constrained to limited regions of 
planet parameter space, it is possible to synthesize results 
from multiple techniques to derive a more accurate and 
complete census o f exoplanets (e.g . iHowar d et alJl2012t 
iGaidos et al.11201'3 : iClanton fc Gaudill2014al lbfh 

Figure Q] shows several existing constraints on the oc¬ 
currence rate of planets around M dwarfs for select 
microlensing, RV, and direct imaging surveys. A di¬ 
rect comparison between these results is complicated 
by the various selection effec ts and observation a l biases 
intrinsic to each survey. In IClanton fc Gaudi ()2J314b|) 
we perform a care ful comparison of the results of the 
iGould et al.1 (120101) microlensing surv ey with the RV 
studies of HARPS dBonfils et al.ll2013D and t he Califor¬ 
nia Planet Survey ICPS: iMontet et alj|20l4) (shown in 
figure [T]) by mapping the observable parameters of the 
population of planets inferred from microlensing to those 
of an analogous population orbiting RV-monitored stars. 
Using the actual detection sensitivities of the HARPS 
and estimated sensitivities for the CPS, we determine 
the expected numbers of planet detections and long-term 
trends for both samples and compare with the reported 
numbers of detections and trends. We find agreement 
with the detections, demonstrating that RV and mi¬ 
crolensing results are consistent, but found that we could 
not explain all the long-term trends. These excess trends 
are either not due to a planetary population, or are due 
to a population of planets to which microlensing is not 
sensitive. 

A similar approach to comparing the results of mi- 
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Figure 1. Planet frequency measurements from several exoplanet discovery surveys of M dwarfs as functions of planet mass, m p , and 
semimajor axis, a. The vertical error bars on the points represent the measurement uncertainty on the planet frequency, while the 
horizontal “error” bars represent the approximate ranges of mass and semimajor axis over which these frequencies are measured. The 
region corresponding to giant planets in the left panel is shaded gray, while the vertical gray line in the right panel represents the 
approximate location of the ice line around a 0.5 Mq star. 


crolensing and direct imaging surveys, where we would 
predict the number of expected planet detections with 
the actual reported number of detections for a given 
imaging survey, would not be as informative as it was 
in the case of microlensing and RV. To zeroth order, the 
“sensitivity function” of direct imaging surveys to detect¬ 
ing planetary companions depends on just two parame¬ 
ters: the planet-to-star contrast ratio and the projected 
separation. Of course, these parameters are, themselves, 
dependent on several system-specific properties — such 
as age, metallicity, spectral type, distance, planetary at¬ 
mospheric composition and structure — that are suscep¬ 
tible to large uncertainties, complicating direct compar¬ 
isons between results of high-contrast imaging surveys 
with those of other detection methods. We define the 
“typical” microlensing planet to be one residing in the 
region of peak sensitivity of microlensing surveys, de¬ 
scribed by the following parameters: a host star mass 
of Mi ~ 0.5 Mq, a mass ratio of q ~ 5 x 10~ 4 , and 
a projected separation of r± ~ Re ~ 2.5 AU, where 
Re ~ 2.5 AU is the typ ical Einstein radius for a 0.5 Mq 
lens (iGould et all 12010 1. This corresponds to a typical 
microlensin g planet mass o f m p = qMi ~ 0.26 Mj up ~ 
Msat- lOuanz et al.l (2012) present the deepest high- 
contrast image for an M star to date, looking at the 
12 — 50 Myr old pre-main-sequence star AP Col, which 
is located 8.4 pc away, for a planetary companion. The 
authors present their 5er detection limits in their figure 2. 
It is clear from this figure that a Saturn-mass planet at 
a projected separation of ~ 2.5 AU (i.e. the typical mi¬ 
crolensing planet) would be undetectable. This suggests 
that most, if not all, of the planets to which microlens¬ 
ing is sensitive will be undetectable to imaging surveys, 
meaning there is probably little overlap between these 
two techniques. 

Thus, we develop new methods in this study to 
constrain the distributions of the properties of wide- 
separation (a > 2 AU) planets by combining the demo¬ 
graphics derived from several individual microlensing, di¬ 


rect imaging, and RV surveys. More specifically, we will 
assume a simple, power-law planet distribution model, 
map the resultant populations into the relevant observ¬ 
ables for each type of survey, and employ the actual de¬ 
tection limits for each survey to determine the expected 
observations. We will then weight the assumed planet 
population by its likelihood associated with each given 
survey, repeating this process many times to derive con¬ 
straints on the parameters that describe our planet dis¬ 
tribution model (i.e. power-law slopes, normalization, 
and the outer cutoff radius of the separation distribution 
function). 

This is the first study to attempt to perform a joint 
analysis combining the results from five different exo¬ 
planet surveys using three independent detection meth¬ 
ods. The constraints derived herein represent the most 
comprehensive understanding of large-separation plan¬ 
ets that orbit the most common stars in our Galaxy and 
will serve as the standard to which models of planet for¬ 
mation and migration must adhere. In § [2j we provide 
brief descriptions of the surveys and their results included 
in our analysis. We perform an order-of-magnitude- 
style analysis in § [3] under several simplifying assump¬ 
tions, such as a single characteristic host mass, distance, 
and kinematic parameters, as well as median sensitivity 
curves to describe each survey. These rough calculations 
allow us to develop the intuition necessary for under¬ 
standing the constraining power of both the individual 
and combined survey results. In §01 relaxing as many of 
these assumptions as possible, we describe the methodol¬ 
ogy we develop to perform a detailed analysis that treats 
uncertainties in the relevant stellar physical and kine¬ 
matic parameters, orbital parameters, planetary evolu¬ 
tionary models, and variations in detection sensitivities 
among, and within, the surveys to which we compare. 
We present and discuss the results of these careful anal¬ 
yses in §0 and follow up with a description of all sources 
of (quantified and unquantified) uncertainty in our anal¬ 
ysis in m We summarize our results and discuss future, 
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related work in § [3 

2. MICROLENSING, DIRECT IMAGING, AND RV SURVEY 
RESULTS 

Here we describe the results from several representa¬ 
tive microlensing, direct imaging, and RV surveys we 
use to constrain the demographics of planets around M 
dwarfs. In this paper, we will make use of the follow¬ 
ing results: 1) the measured fre quency of gia nt planets 
from the microlensing survey of IGould ct all ( 201 0s). 2) 
the slope of the planetary mass- ratio f unction inferred 
from the microlensing survey of iSumi et al.l (120101) . 3) 
the detection of four long-term RV trends measured by 
the CPS TRENDS survey dMontet et al.1120141) . 4) the 
non-detection of any planetary companions in a sam¬ 
ple of 72 young, single M stars im aged by the Planets 
Around Low-Mass Stars (PALMS; iBowler et al.l I2015D 
survey, and 5) the non-detection of any planetary com¬ 
panions around the 16 M st ars imaged by the G emini 
Deep Planet Survey (GDPS: iLafrenibre et alj|2007l) . 

Figure[2]graphically displays three of these constraints. 
Note that the measurement for the CPS/TRENDS sur¬ 
vey dMontet et al.l I2014D differs in this plot relative to 
that in figure [T] Since we are concerned with the long- 
period companions to M dwarfs in this paper, we only use 
the constr aints from the trends that are detected and re¬ 
ported by iMontet et al.l (|2014D (discussed in § 12.31) , and 
thus the data point plotted here is the result only con¬ 
sidering the constraints from their long-term RV drifts. 
Not shown in figure [3 are the two constraint s we con¬ 
sider from the GDPS (lLafreniere et al.l 120071) . as they 
only report a frequency for all stars in their sample (of 
whic h only 16 are M dwarfs; see § 12.21) . and that of the 
ISumi et al.l ( 2010 ) microlensing survey, which is a mea¬ 
surement of the slope (but not the normalization) of the 
planetary mass-ratio function (discussed in the following 
subsection). In the remainder of this section, we describe 
basic sample properties and results of all five surveys we 
consider, and place them into the context of the current 
knowledge of planet occurrences around M dwarfs and 
solar-type stars. 

2.1. Microlensing 

IGould et al.l (|2010f) : The sample of IGould et al.l 
(I2010D is an unbiased selection of 13 high-magnification 
(A > 200) microlensing events. Within this sample, 
there are a total of six planet detections in five systems. 
They find typical values of the planet-to-star mass ra¬ 
tio, Einstein radius, and lens mass of qo ~ 5 x 10 -4 , 
Re,o ~ 3.5 AU(MWM^) 4 / 2 , and M; j0 ~ 0.5 M 0 , respec¬ 
tively. IGould et al.l (l2010l) measure the frequency of ice 
and gas giant planets with mass ratios in the interval 
—4.5 < log q < —2 to be Q = d 2 N p \/(dlogq dlogs) = 
(0.36 ± 0.15) dex -2 at their mean mass-ratio, qo. This 
result is sensitive to a wide range of projected sepa¬ 
rations between s“ 4 x i?E,o < r± < s max i?£,o, where 
Smax ~ (g/10 -4 - 3 ) 1 / 3 , corresponding to de-projected sep¬ 
arations well beyo nd the ic e lines of these systems. In 
their analysis, IGould et al.l ( 2010 ) make the assumption 
that planets are distributed uniformly in log s near the 
Einstein radius and demonstrate that this assumption is 
consistent with the six microlensing detections in their 
sample. 


On its face, this result appears to be in direct con¬ 
flict with those of recent RV surveys, which have mea¬ 
sured frequencies of giant planets around M dwarfs 
that are smaller by more than an order-of-magnitude 
(ICummi ng et alJ l2008l:I.Iohnson et al~ll20 1 0b lllonfils et al.1 

20131: IMontet et, al l 120141 ). Indeed, the publication of 
Gould et al.1 (l2010f) left many uncomfortable with the 
idea that giant planets may be common around M dwarfs 
because the theory of giant planet formation by core ac¬ 
cretion, suppo rted by the RV measu rements, predicts the 
opposite (e.g. lLaughlin et ahl 12004 1). Qualitative argu¬ 
ments invoking the Galactic metallicity gradient or the 
sensitivity of microlensing surveys to longer periods than 
RV surveys seemed likely explanations for this vast differ¬ 
ence in inferred occurrence rates. However, RV surveys 
of M dwarfs have already operated long enough to have 
significant sensitivity, at least in terms of orbital period, 
to the majority of planets detectable by microlensing. 
This motivated a study to statistically compare the con¬ 
straints on exoplanet demographics inferred from these 
two discovery techniques and develop a thorough under¬ 
standing of the conflati on of the various s e lection effects 
intrinsic to them. In iClanton fc Gaudil (|2014al lbt). we 
demonstrate that the giant planet frequencies measured 
by microlensing and RV surveys are actually consistent, 
and without a need to invoke the Galactic metallicity 
gradient. Rather, the steep planetary mass-ratio func¬ 
tion derived from microlensing discoveries implies that 
RV surveys are detecting only the high -mass ta il of the gi¬ 
ant planet population infe rred bv lGould et all ()201(lD . In 
IClanton fe Gaudil (l2014al) . we al so syn t hesize constr aints 
from the m i crolen sing surveys of lGould et al.l (I2010D and 
ISumi et all (I2010D with th ose from the HARPS RV sur- 
vev of lBonfils et all (120131) and determine the frequency 
of Jupiters and super-Jupiters with 1 < m p smi/Mj up < 
13 and periods 1 < P/ days < 10 4 to be fj = 0.0291o q)|, 
a median factor of 4.3 (1.5-14 at 95% confidence) smaller 
than the inferred fre quency of such planets around FGK 
stars of 0.11 ± 0.20 (ICumming et al.lf2008D . Thus, the 
combined microlensing and RV constraints are consistent 
with the generic prediction of core accretion that giant 
planet frequency is correlated with host mass, although 
it remains to be seen if these rates are consistent in detail 
w ith quantitati v e pred icti ons of this model. 

ISumi et al.l (|2010f ): ISumi et~aD (120101) perform a 
likelihood analysis on the 10 microlensing planet dis¬ 
coveries up to their time of publication to determine 
the planetary mass-ratio distribution function, which 
they assume has the form dN p \/dlogq = N 0 q n Q(q — 
®ow)0(<7u P - q), where ® ow = 10 -4 ' 5 and q up = 10~ 2 
are the lower and upper limits on the mass ratio, respec¬ 
tively, and where 0 is the Heaviside step function. The 
authors report a slope of n = —0.68 ± 0.20 but are un¬ 
able to place constraints on the normalization Nq since 
they estimate the relativ e (not absolute) d etection effi¬ 
ciencies as a function of q. ISumi et all (12O1O0 infer a slope 
that is steeper than, but consistent with, that found by 
RV surveys of solar-tv pe stars of Urv = —0.31 ± 0.20 
l| Cummins et al.ll2008l) . 

2.2. Direct Imaging 

Bowler et al. ( [20151 ): The PALMS survey is a direct 
imaging campaign targeting a sample of 122 young M 
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Figure 2. Constraints on planet frequency around M dwarfs for a few of the surveys we consider in this paper as functions of planet 
mass, m p , and semimajor axis, a. The vertical error bars on the points represent the measurement uncertainty on the planet frequency, 
while the horizontal “error” bars represent the approximate ranges of mass and semimajor axis over which these frequencies are measured. 
The region corresponding to giant planets in the left panel is shaded gray, while the vertical gray line in t he right panel repre sents the 
approximate location of the ice line around a 0.5 Mg star. The result plotted for the CPS/ TRENDS survey of lMontet et all ( 2014 1 includes 
only the constraints from the four significant long-term RV trends they detect (see § 12.31 for description). 


stars that have been selected for their youth and proxim¬ 
ity, with the goals of discovering sub-stellar companions 
and determining the frequency of giant planets orbiting 
M dwarfs at separations beyond ~ 10 AU. The median 
spectral type and distance of the full PALMS sample is 
approximately M3.5 and 21 pc, respectively. The median 
age of the sample is 135 Myr and 90% of the targets are 
younger than 620 Myr. Imaging of these 1 22 M dwarfs 
reveal that 44 are close stellar binaries. iBowler et al.l 
(|2015l) do not include these close binaries in their statis¬ 
tical analysis in order to construct a sample that most 
closely matches those of RV surveys, which are routinely 
vetted of such binaries. Their final statistical sample is 
thus composed of 78 single M stars and represents the 
largest imaging search for planets around low-mass stars 
(0.1 — 0.6 Mq) to date. In our analysis, we further refine 
this sample to construct a sample to which we can fairly 
compare the results of other studies (in particular, the 
CPS RV sample of M dwarfs which have been explicitly 
vetted for binaries), eliminating three stars with spectral 
types earlier than M0: TYC 523-573-1 (K7.5), NLTT 
26359 (K5), and TYC 1752-63-1 (K7). We also elimi¬ 
nate three stars which are close visual binaries with an¬ 
gular separations < 1", where there are no stable S-type 
orbits: 1RXS J034231.8+121622 AB, GJ 3629AB, and 
GJ 4185 Aab. This leaves us with a sample of 72 single 
M stars. Some of these stars are actually in wid e bina - 
ries, but following the analysis of IBowler et all ( 2015 ). 
we only examine separations out to which there is a sta¬ 
ble S-type orbit and assume no sensitivity to separations 
b eyond this maximum stable separation. 

IBowler et al.1 (120151) repo rt typical contrasts of 12 — 
14 mag with Keck/NIRC2 (|Wizinow ich et all I200C ) and 
9 — 13 mag w ith Su baru/HiCIAO (IHavano et al.l 20101 : 
iHodapp et all 120081: iSuzuki et al.l 1201 OP at~ an angular 
separation of 1", which correspond to limiting planet 
masses between roughly 0.5 — 10 Mj up at projected sep¬ 
arations between 5 — 33 AU for 85% of their statistical 


sample. Overall, the PALMS survey is most sensitive to 
planets with projected separations between 10 —100 AU, 
demonstrating an ability to detect ~ 5% of 1 Mj up com¬ 
panions, ~ 20% of 2 Mj up companions, and ~ 50% 
of 5 M,T„n compan ions in this range (see figure 24 of 
IBowler et all 120151 for the full PALMS survey sensitiv¬ 
ity map). 

A total of 167 faint point sources (planet candidates) 
are identified aroun d 45 stars from first epoch images. 
IBowler et al.l (12015 ) were able to obtain follow-up imag¬ 
ing and recover 56% of these candidates, showing that 
all are consistent with background stars. The nature 
of the remaining candidates with just a single epoch of 
astrometry is uncertain, however only 8 o f these reside 
at pro jected separations within 100 AU. IBowler et al.1 
(2015f ) do not confirm any planets orbiting the single 
M dwarfs in their sample. These non-detections pro¬ 
vide 95% confidence upper limits on the frequency of 
giant planets (1 — 13 Mj up ) orbiting single M dwarfs be¬ 
tween 10 — 100 AU of < 10.3% (< 16.0%), assuming 
logarithmically-uniform distributions in planet mass and 
semimajor axis and “hot-start” (“cold-sta r t”) planet evo¬ 
lution ary models of lBaraffe et al.l (|2003[ ) (|Fortnev et al.1 
12003). ' 


Bowler et al.l (120151 ) point out that the results of the 


PALMS survey, when compared with similarly large sur¬ 
veys for the direct detection of exoplanets around sam- 
ples of A (e.g. lNiel sen et al. ll2013HVi gan et a l1l2012D and 
FGK stars Ce.g. INielsen fe Closel 120101 : IChauvin et al.1 
I2015D . provide no statistical evidence for a dependency 
of the frequency of long-period, giant planets on host 
star mass. Larger sample sizes in all mass stellar mass 
regimes will be required of future surveys in order to dis¬ 
tinguish small differences in the relative occurrence rates 
of long-period (> 10 AU) giant planets around low- and 
high-mass stars. _ 

Finally, IBowler et al.l (|2015l) report the detection of 
four brown dwarfs, deriving a frequency of brown dwarf 
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companions to single M dwarfs with masses between 
13 — 75 Mj up and separations between 10 — 100 AU of 
2.8]t 1 'g%, consistent with previous estimates for BD com¬ 
panions to M dwarfsiDicterich et al. (1 20 121) an d to solar- 
type primaries iMetchev fe Hillenbrandl (12009(1 . 

In this study, we will employ the sensitivity limits for 
the 72 single M stars targeted by the PALMS survey to 
identify classes of planet populations that are consistent 
with non-detections and compare these with the popula¬ 
tions consistent with the results of the other exoplanet 
surveys from microlen sing, RV, and the Ge mini Deep 
Planet Survey (GDPS; lLafreniere et al.ll2007l) . Th e con - 
trast curves provided in table 5 of iBowler et al.l (|2015f) 
are the deepest contrast curves for each star, but are not 
necessarily the ones employed in their statistical analysis 
(see their § 6.1 for a description of how they select the 
appropriate contrast curves for their analysis). In or¬ 
der to most accurately determine the expected results of 
the PALMS survey for generic (power-law) planet pop¬ 
ulations and compare with the actual reported results, 
we instead use the same set of contrast curves used in 
their statistical analysis (B.P. Bowler, private communi¬ 
cation)._ 

Lafreniere et al. ()2007f) : The GDPS is a direct 
imaging survey of young, nearby stars, with the pri¬ 
mary science driver of constraining the occurrence rates 
of massive, planetary companions at orbital separations 
ranging between 10 — 300 AU. The GDPS sample in¬ 
cludes 85 FGKM stars within 35 pc at ages ranging be¬ 
tween ~ 10 — 6700 Myr. The median values for several 
parameters of the full sample are a distance of 22 pc, 
spectral type of K0, age of ~ 200 Myr, H magnitude 
of 5.75, proper-motion amplitude of 240 mas yr^ 1 , and 
metallicity of [Fe/H] = 0.00 dex (with a sta ndard devi¬ 
ation of 0.21 dex). lLafreniere et al.l (1200711 do not de¬ 
tect any planetary companions, but report a typical 5a 
detection limit sufficient to detect planets with masses 
> 2 Mj up at projected separations between 40 — 200 AU. 
These authors analyze the detection limits of each star 
in their sample, deriving upper limits on the fraction 
of stars hosting at least one planet with m p > 2 Mj up 
to be 0.23 for companions with separations between 
25 — 420 AU and 0.12 for companions with separations 
between 50 — 295 AU at 95% confidence, assuming the 
“hot-start” planet evolutionary models of iBaraffe et alJ 

( 120031 ). ' 

lLafreniere et al.l (|2007l ) provide the 5 a detection limits 
they derive for each star in the GDPS sample in tab¬ 
ular form, which we obtain and use to determine ex¬ 
pected detections from generic planet populations. We 
perform analyses on the subset of M stars in their sam¬ 
ple, which should more closely represent the population 
of microlensing surveys, at least in terms of stellar mass, 
than the more massive stars. There are a total of 16 M 
dwarfs among their full sample. 

While the GDPS is an overall deeper survey than 
PALMS in terms of contrast, the sample sizes differ 
by a factor of ~ 4. Thus, the reported upper limits 
on giant planet frequency at large separations from the 
PALMS survey seem to be more constraining than those 
reported by the GDPS. However, it is difficult to quan¬ 
tify the actual differences in completeness, as the GDPS 
results are quoted over a larger range of semimajor axis 


(25-420 AU) than that from PALMS (10-100 AU). Fur¬ 
thermore, the GDPS result includes the non-detections 
from the higher-mass stars, and it is not clear if the 
planet evolutionary models were incorporated in a consis¬ 
tent fashion between the two surveys. This highlights the 
importance of performing a self-consistent analysis be¬ 
tween surveys using the same detection technique in or¬ 
der to perform a fair comparison of their results. To this 
end, we perform our own analysis of both the GDPS and 
PALMS sample within the same framework, using the 
exact same planet evolutionary models that have been 
interpolated to a given time resolution and extrapolated 
to a larger range of planet masses and to later ages. We 
describe i n detail how we inc orpo rate the evolutionary 
mo dels of lBaraffe et al.l (120031 ) and lFortnev et al.l (120081 ) 
in § 14.41 

2.3. Radial Velocity Drifts 

iMontet et alJ (12011): The CPS T RENDS survey 
dJohnson et al.M2010t IMontet et al.ll2014D targets a sam¬ 
ple of 111 M dwarfs (which they define as having B — V> 
1.44) brighter than V = 13.5, closer than 16 pc, and vet¬ 
ted for known binaries. The sample has a median num¬ 
ber of observations of 29 over a median time baseline 
of 11.8 years, with typical Doppler precisions of a cou¬ 
ple meters per second. With spectral types from M0 to 
M5.5 and masses ranging between 0.10 — 0.64 M 0 (me¬ 
dian of 0.41 M 0 ), the CPS M dwarfs cover a wide range 
in metallicity, —0.81 < [Fe/H] < 0.52, with a slightly 
s ub-Solar median valu e of —0.1 dex. 

IMontet et alJ (|2014D identify four stars within their 
sample that exhibit long-term RV trends, signaling the 
existence of wide-separation companions. In one of the 
cases, curvature is detected in the RV curve, allowing 
IMontet et al.l (I2014D to constrain the companion mass 
(and to a lesser extent, the orbital period) and show the 
l ong-term RV varia tion must be due to a giant planet. 
IMontet et al.l (I2014D obtain images with AO to rule out 
stellar companions and most brown dwarf masses as the 
cause of their observed RV trends. Thus, considering 
the low inf erred frequency of brown dwarf companions to 
M dw arfs (IM etchev fe Hillenbrandll2009l : iDieterich et alJ 
I2012D . IMontet et al.l (|2014D argue these trends are most 
likely due to giant planets rather than lower-mass brown 
dwarfs that are not quite ruled out by their AO obser¬ 
vations. Assuming all their trends are indeed due to 
planetary-mass companions, and including the planets 
with full o rbits tha t have been detected within their sam¬ 
ple, IMontet et alJ (120141) infer a frequency of Jupiters 
and super-Jupiters with 1 < m p /MQ < 13 and or¬ 
bital separations a < 20 AU around M dwarfs to be 
0.065 ± 0.030. This is consistent w ith the constraint we 
derive in lClanton fe Gaudil (l2014bl) from synthesizing mi¬ 
crolensing and RV results of 0.035to;oi9 i n the same mass 
and semimajor axis intervals, and is inconsistent with th e 
result for F GK stars of 0.11 ±0.02 ((Cumming et al.l200§) . 

Although IMontet et alJ (120141 ) infer a frequency of 
Juipters and super-Jupi ters that is consistent with 
IClanton fc Gaudif (|2014bD . it is nevertheless a median 
factor of 1.9 (0.11 — 7.9 at 95% confidence) times larger. 
This is potentially due to a population of very long- 
period super-Jupiters to whic h micr o lensi ng is no t sen¬ 
sitive but is being inferred bv IMontet et al.l (l2014D from 
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these long-term RV trends. It could be that this popula¬ 
tion is responsible for some or all of the sho rt-timescale 
events observed in the microlensing survey of lSumi et al.l 

MJ- 


tions of each of our population parameters (a, /3, A, and 
a ou t) which describe planet populations consistent with 
the results of all surveys. This process is described in 
greater detail in the following sections. 


3. ORDER OF MAGNITUDE EVALUATIONS 


3.1. Individual Constraints 


The primary objective of this study is to constrain 
the properties of wide (a > 2 AU) planetary compan¬ 
ions to M dwarfs by testing if a simple parameterized 
planet distribution model can explain the microlensing, 
RV, and direct imaging data. In order to do this cor¬ 
rectly, we will need to map the parameters of an assumed 
model of the distributions of planet masses and semima¬ 
jor axes into the appropriate observables for each discov¬ 
ery technique and use the actual detection sensitivities 
for each specific survey to determine their expected num¬ 
ber of detections for a given planet population. The full 
procedure will require marginalizing over several stellar, 
Galactic, kinematic, and orbital parameters, but first we 
will forgo several of these steps to estimate — at an order- 
of-magnitude level of precision — the expected results 
for microlensing, direct imaging, and RV surveys. We 
perform this rough calculation first to develop intuition 
on the potential constraining power of each individual 
detection technique, as well as that from their combina¬ 
tion. These evaluations are also a way to validate our 
more careful—but also more complicated and harder to 
interpret—methodology. 

We begin by assuming planet populations that are de¬ 
scribed by four parameters: the slope of the mass distri¬ 
bution function, a, the slope of the separation distribu¬ 
tion function, /?, their normalization, A , and the outer 
cutoff radius of the separation distribution function, a out . 
The populations are thus distributed in planet mass and 
semimajor axis as 


d 2 N pl 

dlogmp dloga 



(1) 


with a maximum semimajor axis of a ou t- We choose to 
normalize our distribution functions at a Saturn mass 
and at a semimajor axis of 2.5 AU, which correspond to 
the mass ratio and projected separation in units of the 
Einstein radius in the region of peak sensitivity for mi¬ 
crolensing surveys, qo = 5 x 10 -4 and So = 1, for a lens 
mass of Mi ~ 0.5 Mq. Note that 2.5 AU is the median 
Einstein radius of an ensemble of simulated microlens¬ 
ing events with parameters consiste nt with those ob¬ 
served for actual samples of events (see lClanton fc Gaudil 
l2014allbl for more information). 

We generate random planet populations, enumerated 
by the index i, each of which is described by a random set 
of parameters: {a*, j3i, Ai, a ou t,i}- We randomly (and in¬ 
dependently) draw values of a and /3 uniformly within 
the range [—4,4], values of log (A/dex -2 ) uniformly 
within the range [—3,3], and values of log (a ou t/AU) uni¬ 
formly within the range [—1,3]. We then map the proper¬ 
ties of these populations into the observables relevant to 
microlensing, RV, and direct imaging surveys, and ap¬ 
ply the appropriate detection criteria to determine the 
expected observations for each survey. Each planet pop¬ 
ulation is then assigned a statistical weight according to 
its likelihood associated with all the surveys to which we 
compare. We thus derive maximum likelihood distribu- 


3.1.1. Microlensing 

For the following order of magnitude comparisons with 
microlensing surveys, we transform between ( q , s) -O- 
(m p , a) assuming a characteristic lens mass of Mi ~ 
0.5 Mq and the median projection angle of a randomly 
oriented, circular orbit where necessary. These assump¬ 
tions will be relaxed in later sections when we perform 
more detailed ca lc ulatio ns. 

iGould et al.l ([2010l) : We analytically compute the 
integrated planet frequency over the mass and semimajor 
axis ranges of 5 Mq < m p < 5 Mj up and 1.4 < a/AU < 
7.2, respectively, for each population. These correspond 
roughly to the ran ges in mass ratio an d projected separa¬ 
tion to which the IGould et al.l ([2010D survey is sensitive 
(—4.5 < logg < —2 and 0.5 < s < 2.5, respectively). 
The integrated frequency is given by the equation 


Qi 



( CL \Pi 

X V 2 5 AU / dlo & m P dlo S a ’ 


ai 


where 


( 2 ) 


_ r 1.4 AU for a ouM < 1.4 AU 

amax ’* \min(7.2 AU,a ouM ) for a ouM > 1.4 AU. 

(3) 

We then weight the parameters of each planet population 
by the likelihoods associated with our calculated values 
of Qi. These likelihoods are assig ned accord ing to the 
likelihood function measured bv IGould et all (2010). to 
which we as sign the variable ££g. 

Similar to I Clanton fc Gaudii (2014b ), we approximate 
Jzfg by first drawing a non-integer number from a Poisson 
distribution wi th a mean va l ue equ al to the number of de¬ 
tections in the IGould et al.l (120101) sample (six) to repre¬ 
sent the average number of planets per star. We then use 
this number to compute the implied pl anet frequency, Q . 
in a similar manner as described in IGould et al.l ( 2010). 
After repeating this process many times, we normalize 
the resultant distribution of planet frequencies such that 
the median value is e qual to the m aximum likelihood 
value reported bv lGould et al.l (I2010D of 0.36 dex -2 . The 
likelihood function we derive from this process is plotted 
in figure [3] and has a median value and 68% uncertainties 
of 0.36^0’J 2 dex -2 . While this is an approximation of the 
actual likelihood function derived bv lGould et al.l (I2010D . 
who report 68% uncertainties of ±0.15 dex -2 , we show 
in § [6] that other sources of uncertainty will dominate 
over that introduced by our approximation of . 

Figure [4] shows the resultant constraints on the param¬ 
eters of populations consistent with this survey. These 
can be understood by considering how changes in the pa¬ 
rameters of our models affect the integrated planet fre¬ 
quency give n by equation ©. We have assumed the mea¬ 
surement of lGould et al.l ([20101 ) uniformly covers a given 
range of q and s, and thus effectively, a given range of 
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Figure 3. Likelihood functions on the measurements made by the various sur veys to which we co mpare. The upper left panel shows our 
approximation to the likelihood function of giant planet frequency measured bv lOould et al.l 112010)), . wh il e the upper right panel shows 

the likelihood function of the slope of the mass-ratio distribution function, JfsiOi measured by Snmi el. al.l (l20IUi). The lower left panel 
shows the like lihood functions of the number of planets existing around stars in both the GDPS HLafreniere"et~ai.1120071) . JS?L 07 i and PALMS 
HBowler et alJ 12015 ). _S?B15; samples, derived from upper limits due to non-detections of planets by these surve ys. T he lower right panel 
plots the likelihood function of the number of RV trend detections, _S?tr, measured by the CPS TRENDS survey IIMontet et al.ll2014l ). The 
median values and 68% confidence intervals are indicated in each plot by a blue, vertical, dotted line and a red, shaded region, respectively, 
and printed at the top of each panel. 


m p and a (although in reality, the range in m p and a 
varies for each system depending on the orbital phase 
and primary mass - this will be addressed in the more 
detailed calculations presented in later sections). The 
black box in figure [5] illustrates this area, over which the 
iGould et al.l (l2010h survey is sensitive. Another impor¬ 
tant detail that must be kept in mind when thinking 
about these constraints is the point at which our planet 
distribution models are normalized. We refer to this as 
the “pivot point” of our models. The black diamond in 
figure [5] shows the location of the pivot point, at a mass 
of Trip = Msat and a semimajor axis of o = 2.5 AU, rela¬ 
tive to the regions of sensitivity for the various surveys. 
Note that the pi vot point lie s i nter ior to the region of 
sensitivity of the IGould et al. (12010?) survey. As we will 
show, this has consequences that are reflected in the con¬ 
straints we derive on the parameters of our planet models 
from comparison with this survey. 

Comparison of ou r planet models with the results of the 
IGould et alJ (I2010D survey essentially does three things: 
1 ) puts a lower bound on a ou t, 2) puts an upper bound 
on A , and 3) constrains a combination of a and A. 
Panels (a), (c), and (e) in figure [4] show that values of 
flout < 1.4 AU are not allowed (although the bin size 
makes it appear as values down to s=s 1 AU are allowed), 
as that would imply there are no planets beyond 1.4 AU, 
which is clearly discrepant with the observations. Pan¬ 


els (a), (b), and (d) demonstrate that there is a maximum 
value of A allowed. This is best understood mathemat¬ 
ically. First consider cases where a ou t > 7.2 AU. At a 
fixed value of A, there are specific values of a and /3 
(that are completely independent of a ou t) that minimize 
the integrals in equation At these minimizing val¬ 
ues, A becomes the only way to adjust G, meaning that 
there is a maximum value of A beyond which the implied 
planet frequency be comes too large to b e consistent with 
the measurement of lGould et al.l ( 2010 ). 

To understand the constraint on the combination of a 
and A , consider what happens to the integrated planet 
frequency at fixed values of a , /3, and a ou t- In this case, 
we can see directly from equation |2]) that Q oc A. Next, 
consider fixed values of a, /?, and A. Larger values of 
flout will result in larger values of G , while smaller a ou t 
result in smaller G ■ However, at fixed /3, a ou t, and A , in¬ 
creasing (decreasing) a will not always lead to increased 
(decreased) values of G ■ This is due to the fact that the 
pivot point l ies interior to the region of sensitivity of the 
IGould et alJ (I2010D survey. For fixed values of G, a 0 ut, 
and A , there are two values of a that will satisfy equa¬ 
tion ([2]). Qualitatively, what is happening in this case, is 
a redistribution of planets from greater numbers of more 
(less) massive planets to greater numbers of less (more) 
massive planets for decreasing (increasing) values of a, 
while conserving the overall planet frequency. This is il- 
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Figure 4. Likelihood contours as a function of pairs of parame¬ 
ters describing planet distribution functions that are found to be 
consistent (in our order of magnitude evaluations) w ith the mea¬ 
surem ent of the integrated planet frequency, Q, by the I Gould et alJ 
(12010 1 microlensing survey. Contours are drawn at levels of 68% 
and 95% of the peak likelihood and are marginalized over all other 
parameters except those being plotted. 

lustrated in panel (b) of figure [4j Taking a slice through 
these contours at a fixed value of A = 0.01 dex -2 results 
in a bimodal distribution in a. If we had a perfect mea¬ 
surement of Q and knew the exact values of /3 and a ou t, 
this distribution would consist of two “spikes,” or delta 
functions, at specific values of a. However, uncertainty in 
the measurement of Q and mutual degeneracies between 
a, /3, and a ou t (see panels c, e, and f) give these modes fi¬ 
nite width. As A increases to a maximum allowed value, 
these modes move closer together and eventually merge 
into a single peak at a value of a that minimizes the con¬ 
tribution of the mass integral in equation ([2]). Panel (d) 
displays a similar behavior in the point constraints of 
P and A , although it is less pronounced than that be¬ 
tween a and A because the range of a covered by the 
iGould et all (120101) survey is smaller relative to that of 
m p . 

Note that even though panels (b), (c), and (f) seem to 
indicate that there is a maximum value of a allowed, this 
is not a real constraint. Instead, it is a consequence of the 
limited range of A we consider, and if we had considered 
smaller values of A , the contours would expand to include 
larger values of a. However, some values of a would not 
be consistent w i th th e distribution of mass ratios seen 
by iGould et al. 1 :201 01 (and similarly, some values of /3 
would not be consistent with their observed distribution 
of projected separations). We address this point in later 
s ections. _ 

iSumi et al.l (120101 ): This survey measures the slope 
of the mass-ratio distribution function, dJV p i/dlogq oc 
q p , reporting a maximum likelihood value and 68% un¬ 
certainties of p = —0.68 ± 0.22. Roughly, the area in 
m p —a space over which this measu rement is made is illus¬ 
trated by the black box in figure [5] ISumi et aD (120101) re¬ 
port their likelihood function of p in their figure 10, which 


Figure 5. Approximate regions of sensitivity for the various sur¬ 
veys we consider. The solid, black l ines b oun d the region for th e 
microlensing surveys of Gould et al. ( 2010) and Sumi et al. ( 2010). 
The gr een region shows the area within which the CPS TRENDS 
survey JMontet et al.|[2014T) has demonstrated an ability to detect 
long-term RV trends, whi le the red regions show the areas within 
which th e GDPS (solid: ILafreniere et al .1 1200711 and the PALMS 
(dashed; [Bowler et al.ll2015l) surveys are sensitive to direct detec¬ 
tions. The solid, black diamond shows the “pivot point” of our 
planet distribution function, i.e. the point at which we choose to 
normalize our distribution function ( m p = Ms a t? & = 2.5 AU). 
Note that this pivot point is lo cated inside the reg ion o f sensitivity 
for the microlensing surveys of Gould et al. (2010) and Sumi et al. 
([201(1 1 and outside the regions of sensitivity for the remaining sur¬ 
veys. This has consequences that are reflected in the constraints 
we derive on the parameters of our planet distribution mo dels for 
the individual surveys (see text, figure [4] and figures [til I ( il l. 

is well described by a Gaussian with mean p p = —0.68 
and standard deviation a p = 0.22. Since q = m p /Mi 
and we are assuming a characteristic lens mass of Mi ~ 
0.5 Mq for this order of magnitude comparison, we have 
that d N n \/d\ogq o c dN„ ]/dloqa, or equivalently, p = a. 
Thus, ISumi et all (|2010|) provide a direct constraint on 
a , and the likelihood function we adopt has the simple 
analytic form JSfsio = exp [— (ct; — p P Y / 2er 2 ]. We show 
a plot of Jz?sio in figure [3] Figure [6] shows the resultant 
constraints on a. In reality, conflation of the planet mass 
and lens mass distribution functions complicate the re¬ 
lationship between a and p , which we consider in later 
sections when we perform more detailed calculations. 

3.1.2. Direct Imaging 

The following order of magnitude comparisons with the 
results of direct imaging surveys use mean detection lim¬ 
its in terms of planet mass and semimajor axis derived 
using the “hot- start” and “co l d-star t” pla net evolution¬ 
ary m odels of iBaraffe et al.l (|2003l) and iFortnev et al.l 
(1200811 . respectively. In more detailed calculations de¬ 
scribed later, we make use of each individual contrast 
curve reported by the surveys, along with the stellar ages 
and distances (and our predicted planet masses) to deter¬ 
mine the expected numbers of detections for direct imag¬ 
ing surveys, again using both “hot-” and “cold-start” 
models. 
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Figure 6. Likelihood contours as a function of pairs of parame¬ 
ters describing planet distribution functions that are found to be 
consistent (in our order of magnitude evaluations) with the mea- 
sur ement of the slope of the mass-ratio distribution function by 
the lSumi et al.l (120101) microlensing survey. Contours are drawn at 
levels of 68% and 95% of the peak likelihood and are marginalized 
over all other parameters except those being plotted. 


models. Specifically, we fit (by eye) the contour levels, 
which correspond to the detection probability of a given 
planet mass, at a constant semimajor axis (50 AU) as a 
function of log m p , finding 


U d et(log rrip) 


( 0.7 log (TOp/Mj U p) + 0.07 for “hot-start,” 
(0.5 log (mp/Mjup) + 0.01 for “cold-start,” 

( 5 ) 


for planet masses 1 < m p /Mj up < 13, and assume 
that this probability is independent of a. We set 
Pdet(\ogm p ) = 0 for m p < Mj up (since iBowler et all 
120151 measure the frequency of Jupiter- and super- 
Jupiter-mass planets) and further assume this form of 
Pdet(log m p ) is the same for each of the 72 single M stars 
we consider from the PALMS sample. 

Under these assumptions, for a given planet popula¬ 
tion, the number of expected direct imaging detections 
per star, Nbisj, for the PALMS survey is thus deter¬ 
mined analytically as 


Nb15,i =A 


/•loga„ 


(•log (13 Mjup) 


/log (10 AU) Jlog(l Mj U p) 


m p 

M Sa 


x ( 2 5 °au ) p Aet{\ogm p )d\ogm p dloga , 


( 6 ) 


where we now have 


IBowler et al.l (2015 ): We compute the expected 
number of planet detections per star for the PALMS sur¬ 
vey for a given planet population as 


/'log Umax, 

=Ai / 

i Aog TYlp t max 

log flmin 

J log ITlp,min 



X ( 2.5 AU. 

J ‘Wet (log r 


A/s a 


x d log m p d log a 


( 4 ) 


where a max ,i = max (a m i n , a ou t,i) and where the function 
$det(logm p , log a) represents the detection limits of the 
PALMS survey in terms of logm p and logo. The inte¬ 
gration limits of equation (3]) are set by the region over 
which we assume our planet distribution function holds. 
We adopt the values a m i n = 1 AU, m Pjm i n = 1 M ffi , and 

max = 13 /Ujup- 

We assume a very simple form for the de¬ 
tection sensitivity function of the PALMS sur¬ 
vey of $det(log m pi logo) = P de t(logTOp)0(loga - 
log U m i n; det,img)0(log Umax,det,img log Cl), where 0 is the 
Heaviside step function, Pdet(logm p ) is the probability 
of detecting a planet of a given mass (which we as¬ 
sume to be independent of a in the range of semimajor 
axes over which we integrate) and where a m i n ,det,img and 
Umax,det,img are the minimum- and maximum-detectable 
semimajor axes, respectively, for the survey. Since the 
main result of the PALMS survey to which we will com¬ 
pare this order-of-magnitude estimate is the upper limit 
on the occurrence rate of giant planets with semimajor 
axes between 10 — 100 AU, we set a m i ni det,img = 10 AU 

and Umax,det,img = 100 AU. 

We determ i ne Pd ^tilogrrir,) by examining Figure 24 of 
IBowler et al.l (I2015T) . which displays the mean survey sen¬ 
sitivities for both hot- and cold-start planet evolutionary 


_ Jmax(10 AU, a ou t,i) for a ou t,i < 100 AU, 

Omax.i - 1 100 Au for a ouM > 100 AU, 

( 7 ) 

and where Pdet(logm p ) is given by equation (O. The 
total number of expected detections for the full sample 
is then Wot, bis ,i = 72/Vbi5,j. We assign a statistical 
weight to each planet population according to the like¬ 
lihood function we define to be jS?bi 5 = P(0|Wot,Bi5,i)> 
where P(0|Wot,Bi5,i_) = exp(-Wot,Bi5,i) is the Poisson 
probability of detecting zero planets among the 72 single 
M stars in the PALMS sample if the true number of de¬ 
tectable planets around such stars is actually Wot,Bi5,i- 
We plot this likelihood function in figure [3l where the 
abscissa in this case corresponds to iVtot,Bi5,i- 

Figure 0 shows the resultant constraints on the pa¬ 
rameters of populations we find to be consistent with 
this survey. As can be seen in equation 3]), there are 
degeneracies between a, /3, and A which result in no 
constraining power on combinations of these parameters 
from direct imaging results. The only real constraints 
(although weak) allowed by the PALMS survey individ¬ 
ually are on a ou t- Panels (a), (c), and (e) demonstrate 
that there is a maximum allowed value of a ou t for a given 
value of a, /?, and A. The shapes of these contour lines 
are straightforward to understand. For a ou t < 10 AU, 
which is the minimum separation at which PALMS is 
sensitive to planets, there is no constraint on a ou t because 
all planets, regardless of their distribution (and overall 
occurrence rate) in mass and semimajor axis, would be 
undetectable by PALMS and thus consistent with the 
reported non-detections. 

Panel (a) illustrates that for a small enough normal¬ 
ization, aout can be arbitrarily large (within the range of 
Uout we consider). Similarly, there are values of a and 
/3 below which there is no constraint on a ou t) where the 
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Constraints: 


68% Contours 
95% Contours 


Bowler+ (2015) 


- Hot Start 

— Cold Start 




log(A [dex 2 ]) og(a nu , [AU]) a 


Figure 7. Likelihood contours as a function of pairs of parameters 
describing planet distribution functions that are found to be consis¬ 
tent (in our order of magnitude evaluations) with the non-detection 
of a ny planeta ry co mpanions in the PALMS direct imaging survey 
llBowler et al.ll2015fl . Contours are drawn at levels of 68% and 95% 
of the peak likelihood and are marginalized over all other param¬ 
eters except those being plotted. Note that the contours enclose 
the bottom left regions of the panels. 

planets are distributed in such a way that they are un¬ 
detectable. Since in the above-described order of magni¬ 
tude calculation, we have assumed the sensitivity of the 
PALMS survey is the same for each star in their sample, 
the transition from having a constraint on a ou t to hav¬ 
ing no constraint is quite sharp. As we will demonstrate, 
this will not be the case when we perform more careful 
calculations that make use of e ach individual sensitivity 
curve reported bv lBowlcr et al.1 (12015H . Variations in the 
limiting contrast (and thus, planet mass) at a given sepa¬ 
ration among the stars in the PALMS survey will provide 
more sensitivity to these parameters, particularly a, as 
these detection limits will sample our planet distribution 
functions at a range of masses and semimajor axes that 
is not probed by the mean sensitivity curve we employ 
for this order of magnitude comparison. 

Finally, we note that there is not much difference be¬ 
tween the constraints obtained from the “hot-start” mod¬ 
els versus those obtained using the “col d-start” models. 
We explain the reasons for t his in § 15.21 

Lafreniere et al. (|2007fl : Our comparison with the 
GDPS is similar to that described above for the PALMS 
survey, where the number of expected planet detec¬ 
tions per star for a given planet population is given 
by equation d4]). Here, we assume a different form for 
$det(logm p ,logo) to represent the average detection ef¬ 
ficiency of the GDPS in terms o f logw „ and logo. In 
their figure 4, lLafreniere et al.1 (12007T ) present a me¬ 
dian detection limit for a star at the median distance 
(~ 22 pc) of their sample. We assume the relation 
4>det(log m p , logo) = 0(logTOp-logm p , m in,LO7)0(loga- 
loga min ,LO7)0(loga maXi LO7 - log a), where m p , m in , L 07, 
Umin,L07: and «max,L07 are the minimum-detectable 
planet mass, and the minimum- and maximum- 
detectable semimajor axes, respectively, for the survey. 


There is a maximum-detectable semimajor axis that is 
set by the field of view (FOV) of the GDPS observations 
of ~ 11" in radius, which corresponds a projected sepa¬ 
ration in physical units of roughly rj_ ~ 250 AU at their 
median target distance, and thus we adopt a max ,L 07 ~ 
290, implicitly assuming a randomly-oriented, circular 
orbit. 

We infer values of m p , m i n .L 07 by taking the max¬ 
imum contrast to which the GDPS is sensitive of 
Amag = 16.3 in the CH4-short filter (given by the 
representative median sensitivity curve plotted in their 
figure 4) and converting this to a planet mass us¬ 
ing their median distance, an age of 100 Myr, and a 
host absolute H-band magnitude of 4.0 (assumed to be 
equ ivalent to its magnitu de in the CH4-short filter— 
see lLafreniere et all 120071 ). along with the predictions 
of planet evol utionary model s. For the “hot-start” 
models of iBaraffe et all (12003K we derive a minimum- 
detectable mass of mn.mir, j, 0 7 .hot ~ 1.4 M.Tnn, and for 
the ” cold-start” models of iFbrtnev et al.l (1200811 we de¬ 
rive a minimum-detectable mass of m Piin i n ,L07,coid ~ 
3.2 Mj U p. We adopt a minimum-detectable semimajor 
axis of a m ; ni L 07 ~ 12 AU (which corresponds to a me¬ 
dian value of rj_ jm i n> L 07 ~ 10 AU for a ra ndomly-oriented, 
circula r orbit) directly from figure 4 of [Lafrenier e et ahl 
(|2007f) . We again make the assumption that the resul¬ 
tant form of 4>det(logm p , logo) is representative of all 16 
M dwarfs in the GDPS sample. 

Under these assumptions, equation (j4j) takes the form 


^07,1 =Ai 

/•log (ctmax,i) 

/■log (13 Mj up ) / 

m p \ ai 

/log (12 AU) 

-/log Trip, min,L07 \ 

M Sat ) 


X I 

f a \A 

dlogm p dloga 


^2.5 AU/ 

1 


where we have 


^max 


,i 


and 


f max(12 AU, a ont ,i) for a ou t,i < 290 AU, 
(290 AU for a out * > 290 AU, 

(9) 


f 1.4 Mj up for “hot-start”, 
(3.2 Mj up for “cold-start”. 


( 10 ) 


For a given planet population, the number of expected 
planet detections per star, IVL 07 ,i, for the GDPS M stars 
can then be determined analytically via equation (l8l) . 
The total number of detections for the 16 M stars in the 
GDPS is simply Wot,L07,i = 16At L07 ,i. As before, we 
assign a statistical weight to each planet population ac¬ 
cording to the likelihood function we define to be 2 z?lo 7 = 
U(0|Atot,L07,iX where P(0|Atot,L07,i) = exp(-A t ot,L07,i) 
is the Poisson probability of detecting zero planets among 
the 16 M stars in the GDPS sample if the true num¬ 
ber of detectable planets around such stars is actually 
Atot,L 07 ,i- We plot this likelihood function in figure 0 
where the abscissa in this case corresponds to -/Vtot,L07,i- 
Figured] shows the resultant constraints on the param¬ 
eters of populations we find to be consistent with the 
GDPS. The interpretation of these constraints (and non¬ 
constraints) is identical to that described for the PALMS 
survey above. Note that although the GDPS is, in gen¬ 
eral, deeper than the PALMS survey, their constraints 
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Figure 8. Likelihood contours as a function of pairs of parameters 
describing planet distribution functions that are found to be consis¬ 
tent (in our order of magnitude evaluations) with the non-detectio n 
of any planetary companions in the GDPS ijLafreniere et afll2007D . 
Contours are drawn at levels of 68% and 95% of the peak likelihood 
and are marginalized over all other parameters except those being 
plotted. Note that the contours enclose the bottom left regions of 
the panels. 


are somewhat weaker due to the fact that the PALMS 
sample size is larger than that of the GDPS by a factor 


of 4.5. 


3.1.3. Radial Velocities 

iMontet et al.1 (l2014f) report the measurement of four 
significant RV trends in the CPS M dwarf sample 
due (most probably) to planetary companions, demon¬ 
strating an ability to measure accelerations of v ~ 
lms -1 yr^ 1 over long time baselines (~ 10 yr). In 
order to estimate the number of long-term RV trends, 
expected for the CPS M dwarfs for a given planet 
population, we use the m edian sensitivity we derive in 
IClanton fe Gaudil (l2014bll for the CPS sample to deter¬ 
mine the minimum-detectable acceleration as a function 
of semimajor axis, ^thr(a)- Mathematically, the num¬ 
ber of expected detectable trends is quantified by equa¬ 
tion @ with -/Vdet.img.i -t N tl> i and <L det (log m p , log a) -t 
0(P — T)Q(v — ilthr), where T is the time baseline of ob¬ 
servations for the CPS TRENDS survey. This integral 
cannot be evaluated analytically, so we instead evaluate 
it numerically with some assumptions. 

First, we assume all planets described by our distribu¬ 
tion function are on randomly oriented, circular orbits. 
In order to compute the magnitude of the RV trends 
these planets would induce, v, we employ Kepler’s third 
law (in the limit M; m p ), which yields the period 

< n > 

For our typical parameter values, the pe¬ 
riod described by equation m evaluates to 
P « 5.6 yr(M;/0.5 M Q )~ 1 / 2 (a/ 2.5 AU) 3/2 . 


Equation m in conjunction with the standard radial 
velocity equation yield the velocity semi-amplitude for a 
circular orbit as: 


K(Mi,a,m p ,i) 


G \ 1/2 

—j m p s rnt, 


( 12 ) 


where i is the inclination of the orbit. We adopt 
a convention wherein the inclination angle is mea¬ 
sured relative to the plane of the sky, such that 
such that i = 0° corresponds to an orbit viewed 
face-on, while i = 90° corresponds to an orbit 

viewed edge-on. Equation m evaluates to K m 
7.6 m s _1 (M^/0.5 M 0 ) _1 / 2 (a/2.5 AU) _1 / 2 (m p /Mg at ) sini 
for our typical parameter values. 

The magnitude of the acceleration of a host star due to 
a companion on a circular orbit is v = (2nK cos Mq)/P 
(where Mq is the mean anomaly, which describes the 
orbital phase at a given epoch), which combined with 
equations (HHH21) takes the form 


v(a,m p ,i,M 0 ) 



sin i cosM 0 . 


(13) 


Typically, this corresponds to an acceleration of 
v w 8.6 ms -1 yr _1 (a/2.5 AU) _2 (m p /Ms a t) sinicosM 0 . 
When computing these accelerations, rather than 
marginalizing over unknown orbital parameters, we 
will assume that the inclination and mean anomaly of 
each system are at their median respective values of 
sinz| med ~ 0.867 and cosM 0 | med = 0.5 for the purposes 
of this order of magnitude calculation. 

To determine which planets produce accelerations 
which are actually detectable by the CPS, we estimate 
the minimum-detectable acceleration as a function of 
their orbital per iod (and thus, imp li city, a) in the follow¬ 
ing manner. In IClanton fe Gaudil (l2014bl) . we estimate 
the sensitivity of the CPS to each of their M dwarfs in 
terms of a signal-to-noise (S/N) ratio. We adopt these 
same sensitivities for this order of magnitud e evaluation 
to infer v t hr (P) ■ In IClanton fc Gaudil (I2014alf , we also de¬ 
rive an expression for the S/N at which an RV survey can 
detect planetary companions as a function of the mag¬ 
nitudes of the signal K and measurement uncertainties 
cr, the time baseline of observations T, and total number 
of observations N. We found the “phase-averaged” S/N, 
which we designate as Q, to be 



(14) 


for which we fin d the CPS sample to have a median value 
°f Qmcd = 8.26 ([Clanton &: Gaudill2014bll . 

Solving equation m for K yields the minimum- 
detectable velocity semi-amplitude, K m ; n , for a given 
value of Q , which we then use to determine the minimum- 
detectable acceleration as i>thr(P) = (27rAT m i n )/P. We 
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derive this to be 



which is implicity a function of a via Kepler’s third law 
(equation 1111) . Since we are aiming for a rough estimate, 
we take the median parameters of the CPS M dwarfs to 
determine u t h r and assume the resulting threshold as a 
function of a is the same for each of the HIM dwarfs in 
their sample. The median parameters for the CPS sam¬ 
ple are iV m ed = 29, er me d = 4.5 m s -1 (indues instrumen¬ 
tal noise and stellar jitter), and T me d = 11.8 yr. We plot 
this minimum-detectable acceleration as a function of log 
semimajor axis in figure [9] The dotted blue, vertical line 
marks the log semimajor axis that roughly corresponds 
to the median time baseline of observations, T me d, for 
the CPS M dwarfs, log ay ~ 0.59, assuming the median 
stellar mass of the CPS M dwarfs of M* ~ 0.41 Mq. As 
a increases beyond ay, Uthr becomes an extremely weak 
function of P (and thus a). Expanding the Maclaurin 
series of sin 2 (x) = x 2 — x 4 /3 + 0(x 6 ), where x = nrT/P, 
and substituting this expression into equation (fl5j) . we 
find 



in the limit of T -C P. This evaluates to u t hr ~ 
2.9 m s^ 1 yr^ 1 at the median parameter values for the 
CPS M dwarfs listed above. We plot this in figure [9] as a 
red, dashed, horizontal line, verifying that it is a good ap¬ 
proximation for equation usd, shown by the solid, black 
curve, for a^> ar- 

For each planet population, we numerically evalu¬ 
ate equation making use of equations CD. CD. 
and GSD. and adopting $det(log m p , log a) —> 0(P — 
T)0{v Uthr), Umin — 0.1 AU, fl ma x,z — max(a m i n , n 0 ut,z)> 
nip,min = 1 M®, and m p , max = 13 M Jup . This yields 
the expected number of long-term RV trend detections 
per star, N tr .i , for the CPS TRENDS survey, assum¬ 
ing their median sample parameters. The total num¬ 
ber of trends we expect the CPS to detect for a given 
planet po pulation is then Wot,try = Ill-Wry- Since 
iMontet et al.l ( 2014 (1 report the detection of four signifi¬ 
cant long-term RV trends, we adopt a likelihood function 
of the form Jzf tr = -P (4 1 Wot,try), where P(4|Wot,try) = 
(Wot,try) 4 exp(— Wot,try)/4! is thc Poisson probability of 
detecting exactly four trends from the 111 M stars in 
the CPS sample if Wot,try is the true mean number of 
detectable trends within the sample. We plot this likeli¬ 
hood function in figure 0 where the abscissa in this case 
corresponds to Wot,try- 

Figure |TU] shows the resultant constraints on the pa¬ 
rameters of populat ions consistent with the results of 
IMontet et al.l (2014). There is not much orthogonal¬ 
ity between t hese c onstraints and those supplied by the 
iGould et al.l (120101) microlensing survey; the RV trends 
also put a lower bound on a ou t and constrain a combi¬ 
nation of a and A, but do not put an upper bound on 
A. The lower bound on a ou t (panels a, c, and e) comes 



0.0 0.5 1.0 1.5 2.0 2.5 3.0 
log(a [AU]) 


Figure 9. The minimum-detectable acceleration as a function of 
semimajor axis at a S/N of Q = 8.26 for an RV-monitored target 
with N = 29 observations over a time baseline of T = 11.8 yr, 
measurement uncertainties of a = 4.5 m s 1 , and a stellar mass 
of M* = 0.41 Mq. These values are the median parameters of 
the CPS M dwarf sample. The dotted blue, vertical line marks 
the log semimajor axis that roughly corresponds to T = 11.8 yr, 
log ~ 0.59, assuming M* = 0.41 Mq. In the limit a 'JP ay, 
i’thr becomes a very weak function of the orbital period (and thus 
a), and is approximately 2.9 m s —1 yr -1 for the chosen parameter 
values (shown by the red, dashed line). 

from the fact that we require the period of planets that 
we expect to be detected as a trend to be longer than 
the time baseline of observations, i.e. P > T. Since we 
use only a median sensitivity to be representative of the 
full CPS M dwarf sample, this lower bound on a ou t is 
sharp at a = ay = 3.85 AU. Later in this paper, we 
will make use of each individual sensitivity limit and the 
lower bound placed on a ou t will not be as sharp, as the 
stars in the CPS sample have been monitored for a range 
of T. The observed anti-correlation between a and A in 
panel (b) is a reflection of the fact that the pivot point of 
our planet distribution function is outside the region of 
sensitivity for the CPS trends, as illustrated in figure [5] 
At fixed aout and /?, increases (decreases) in A must be 
compensated by decreases (increases) in a to preserve 
the number of observed trends. 

The apparent constraints on the combination of a and 
/3 in panel (f) are not actual constraints; if we were to 
consider a larger range of A, this full space would be 
allowed. Similarly, there would not appear to be a max¬ 
imum allowed value of cn, as implied by panel (c). 

3.2. Combined Order of Magnitude Constraints 

Individually, the results of any single exoplanet dis¬ 
covery survey are unable to place precise constraints on 
each of the four parameters of our assumed planet distri¬ 
bution function. The results of single surveys are limited 
to small regions of planet parameter space (see figure 0, 
but we can obtain real constraining power by simultane¬ 
ously comparing with all surveys. 

As s hown in the preceeding sections, the IGould et alJ 
(|2010l) microlensing survey constrains the combination of 
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Figure 10. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent (in our order of magnitude evaluations) with the detec¬ 
tion of four significant, long -term RV trends in the CPS M dwarf 
sample JMontet et al.l[20141) . Contours are drawn at levels of 68% 
and 95% of the peak likelihood and are marginalized over all other 
parameters except those being plotted. 


Figure 11. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent (in our order of ma gnitude evaluations) wit h the results 
of the microlensing surveys of I Gould et al. ( 2010) and Sum i et al.l 
(2010). Contours are drawn at levels of 68% and 95% of the peak 
likelihood and are marginalized over all other parameters except 
those being plotted. 


a and A. Since the iSumi ct~ahl 12010 *) survey basically 
fixes a, th e combination of th is constraint with the re¬ 
sults of the lGould et al.l (120101) survey then fixes A. This 
is illustrated in panel (b) of figure [Tl] which displays the 
combined constraints on the para meters of our planet 
dist ribution fun c tion f rom both the lOould et al.l (|2010f ) 
and ISumi et all (|2010l ) microlensing surveys. However, 
we find that a degeneracy between /3 and A remains 
(panel d), and there is still no constraining power on 
Gout; other than the previously described lower limit of 
aout 1-4 AU. 

With such tight constraints on a and A, t he re¬ 
sults of direct ima ging surveys (lLafreniere et al.l 120071 : 
iBowler et al .1120 15 ) and those of the CPS TRENDS sur- 
vev dMontet et aLH2014i) become sensitive to the values 
of /3 and a out in complementary ways. As we will demon¬ 
strate, the RV trends exclude small values of /3 and put 
a lower bound on a ou t, while direct imaging surveys tend 
to exclude large values of /3 and limit large values of a ou t 
(to a certain point that is dictated by the FOV of the 
direct imaging surveys). 

Figure [12] shows our constraints from the combination 
of microlensing and direct imaging surveys. Panels (d), 
(e), and (f) demonstrate that inclusion of direct imaging 
results excludes large values of /3 (for a ou t ^ 10 AU where 
there is no sensitivity from direct imaging) which would 
predict (at the essentially fixed values of a and A) large 
numbers of direct detections that are inconsistent with 
the observations. Panel (e) shows that as /? decreases, 
aout must decrease, up to the maximum separation at 
which these surveys are able to detect companions (~ 
300 AU for the GDPS), set by their FOV. 

Figure [T3l shows our constraints from the combination 
of results from microlensing surveys and the trends de¬ 
tected by the CPS TRENDS survey. Panels (d), (e), 


and (f) demonstrate that inclusion of RV trend detec¬ 
tions excludes small values of /3 that would predict (at 
the essentially fixed values of a and A) too few long-term 
tr end detections to be consistent with the observations 
of iMontet et al.1 (120141) . Panels (a), (c), and (e) show 
that there is a lower bound on a ou t, interior to which we 
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Figure 12. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent (in ou r order of magnitude evaluations) with the results 
of microlensing (Gould ct al. 2010; Sumi et al. 2010) and direct 
imaging (Lafrcnierc et al. 2007; Bowler et al. 2015) surveys. Con¬ 
tours are drawn at levels of 68% and 95% of the peak likelihood 
and are marginalized over all other parameters except those being 
plotted. 
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would predict no planets to be detected as trends. Since 
the acceleration of a host star induced by the presence of 
a planetary companion depends strongly on orbital sep¬ 
aration (v oc a -2 ), the sensitivity of RV surveys to trend 
detections rapidly declines with a. Thus, there is no con¬ 
straining power on the upper limit fora ou t by inclusion 
of the RV results of I Mont et et al.l (|2014f) . 

Finally, combining the results of all microlensing, di¬ 
rect imaging, and RV surveys we consider, we derive 
the constraints displayed in figures [TJ] (“hot-start mod¬ 
els”) and m ( “cold-start models”). For reference, we 
have also included the constraints obtained by combin¬ 
ing the above-described combinations of survey results. 
These constraints are derived from order of magnitude- 
level comparisons, but suggest that our simple, power- 
law planet distribution function can simultaneously ex¬ 
plain the results of multiple exoplanet surveys, warrant¬ 
ing more detailed and careful calculations (described in 

§[3J. 

We plot the corresponding one-dimensional likelihood 
functions on each parameter, marginalized over the other 
three parameters, in figure 1161 and summarize the me¬ 
dian values (and their 68% confidence intervals) in ta¬ 
ble |T| for both “hot-” and “cold-start” planet evolution¬ 
ary models. 

We have shown that our order of magnitude compar¬ 
isons reveal that there is overlap in parameter space from 
the results of microlensing, RV, and direct imaging sur¬ 
veys. We are also interested in determining the probabil¬ 
ity of the resultant parameter distributions given by each 
survey individually. That is, how much does the likeli¬ 
hood function from a given survey vary over the 68% 
confidence interval determined from the combination of 
all surveys? In figures [T7] (“hot-start”) and [T8] (“cold- 
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Figure 13. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent (in our order of magnitude evaluations) with the re sults 
of microlensing surveys (Gould et al. 2010; Sumi et al. 2010) and 
the RV trends detected by the CPS TRENDS survey ( Montet et al.| 
120141) . Contours are drawn at levels of 68% and 95% of the peak 
likelihood and are marginalized over all other parameters except 
those being plotted. 


start”), we plot the normalized likelihood functions (i.e. 
the likelihoods normalized by the maximum likelihood 
value for each individual survey) of given parameters 
across their 68% confidence intervals, assuming the re¬ 
maining three parameters are at their median values (see 
tabic [U). These figures show that the likelihoods of these 
parameters from each individual survey are consistent to 
within ler (i.e. their normalized likelihood is greater than 
~ 0.32) across much of their inferred 68% confidence in¬ 
tervals, and with normalized likelihoods > 0.85 at the 
median value for all surveys (indicated by the vertical, 
grey lines). This demonstrates that this overlap in pa¬ 
rameter space we identify is not improbable, but has a 
rather good chance of being real. In other words, these 
order of magnitude evaluations indicate that the results 
of microlensing, RV, and direct imaging surveys can be 
plausibly explained by a single planet population mod¬ 
eled with power-law distributions in mass and semimajor 
axis. 

Note that while the direct imaging results, particu¬ 
larly those of the PALMS survey, in the “cold-start” case 
(figure 1181) seem to be more constraining than those in 
the “hot-start” case (figure H71). this is not actually the 
case. The median value of a out in the “cold-start” case is 
larger than the minimum-detectable semimajor axis for 
the PALMS survey (10 AU), while this is not so for the 
“hot-start” case (refer to table [T]). This is not surpris¬ 
ing, as one would expect a larger median value on a out 
for the “cold-start” models relative to the “hot-start” 
models. However, it does make the interpretation of fig¬ 
ures m and HHa bit complicated. Since the values of 
the parameters over which we marginalize in each panel 
of these figures are set to the their respective median 
values, direct imaging surveys are fully consistent (nor¬ 
malized likelihood of one) over the full 68% confidence 
interval of the parameter plotted in each panel, as zero 
detections are predicted for a ou t < 10 AU, with the ex¬ 
ception of the top-left panel displaying the likelihood of 
aout itself. Thus, since the median value of a ou t is larger 
than 10 AU for the “cold-start” models, the number of 
expected direct imaging detections is non-zero, giving 
rise to the constraints within the 68% confidence inter¬ 
vals each parameter displayed in figure [TH 

4. METHODOLOGY 

The order of magnitude-style calculations presented in 
the last section are useful for developing intuition about 
the ways in which individual surveys and their combina¬ 
tion can constrain the distribution of properties of plan¬ 
ets. In order to obtain more robust constraints, we now 
relax as many assumptions from those calculations as 
possible. In this section, we account for variations in 
stellar properties (mass, distance, age), Galactic prop¬ 
erties (velocity dispersion, stellar number density along 
a given sight line), planetary orbital parameters, planet 
evolutionary models (“hot-start”, “cold-start”), and de¬ 
tection sensitivities for each star/target in each survey. 

The general procedure we adopt is as follows. We begin 
by generating a population of planets from our simple, 
power-law distribution function, assigning randomly ori¬ 
ented orbits and host star properties. Next, we derive the 
expected observables for each planet in the population 
for microlensing, RV, and direct imaging observations. 
Applying the sensitivity limits of the various surveys, we 
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Figure 14. Likelihood contours as a function of pairs of parameters describing p lanet dist ributi on f uncti ons that are found to be consistent 
(in o ur order of magnitude eval uations) with the r esults of microl ensing (G10 -Gould et al. 2010: SlO-Sumi et al. 2010ft . direct imaging 
fL07- [Lafreniere et al. 1(20071 : B15 4Bowler et al .11 20 151), and RV (M144Montet et alJl2014l) surveys. Here, we have utilized the predictions of 
the “hot-start” planet evolutionary models of[Baraffc et al. ( 20011). Contours are drawn at levels of 68% and 95% of the peak likelihood 
and are marginalized over all other parameters except those being plotted. 


Table 1 

Median values and 68% uncertainties on t he parameters of ou r pla net distribut i on mo dels for our order of magnitude com parisons with 
the results of the microlensing surveys bv IGould et al.l (2010) and Sinni ~ all (120 1 Oi l, the GDPS (Lafrcnicrc ct al. 200® and PALMS 
HBowler et al .if^OTETl direct imaging surveys, and the CPS TRENDS survey of long-term RV trends IIMontet et al.M2014fl . 


Planet Evolutionary 

Models 

Median 

a 

Values and 

p 

68% Uncerta 

A [dex -2 ] 

linties 

&out [AU] 

“Hot-Start” 

(Baraffe et al. 2003) 


0.49±?;11 

0 26”*” 0 ' 21 
u - zo —0.18 

Q 7+ 62 

' —3.6 

^^Tjold- StarP^™ 

(Fortnev et al. 2008) 

-0.58i°;?° 

o-3i±?;|| 

r» 9(7+0.20 
u ' z ‘-0.18 

VI+105 

X1 -5 


predict the results of each survey given the specific planet 
population under consideration, comparing these predic¬ 
tions with the actual results to assign a likelihood to 
the population. This process is repeated many different 
planet populations until we have adequately sampled all 
the relevant parameter space. We describe this method¬ 


ology in detail in the remainder of this section. 

4.1. Generating Planet Populations 
4.1.1. Masses and Separations 

In a given realization, *, corresponding to a single 
planet population described by {crj, /3j, + , a ou t,+ we 
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Figure 15. Identical to figure fill but utilizing the predictions of the 

generate an ensemble of N ens = 10 6 planets enumerated 
by the index j. Each planet is described by a mass, 
and a semimajor axis, a,j, drawn from the distri¬ 
bution function 


d 2 N p \ 

dlogrrip dloga 




(17) 


The average number of planets per star in each realiza¬ 
tion within the range of planetary masses and separations 
we consider is then given by 


/■log (a„ut,i) /-log (m p , ma *) d 2 N } 

iiog(a in ) iiog(m p , min ) d log m p dloga i 
x d\ogm p dloga . (18) 


For all realizations, we fix the minimum semimajor axis 
at Oi n = 0.1 AU, and the minimum and maximum planet 
masses are set to m Pim i n = 1 M® and m Pimax = 13 Mj up . 


“cold-start” planet evolutionary models of IFortnev et all ; 62008 1). 

Mapping the masses and semimajor axes of planets 
belonging to such populations into the observables for 
various types of surveys requires generating orbits. We 
assign a random orbit to each planet by drawing values 
from prior distributions of eccentricity, e,; >7 , inclination, 
ii j , argument of pericenter, Ui j, and mean anomaly, 
Note that we do not need to assign a longitude of 
ascending node, to our planets since the projected 
separation (required for microlensing and direct imaging) 
is independent of this angle, as are the relevant RV ob¬ 
servables. Ecc entricities are drawn fro m the distribution 
determined bv iZakamska et abl (120111) . which consists of 
two components. The first, representing 38% of the pop¬ 
ulation, is characterized by circular orbits, while the sec¬ 
ond, repres enting “dynamical l y acti ve” planets, follows 
the form of lJuric fe Trcmaind (|2008f ) given by 


dN pl 

-oc e exp 




4.1.2. Orbits 


de 


(19) 
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Figure 16. Relative likelihoods fo r each parameter, mar ginalized over the other three paramet ers, derived from our o rder of magnitude 
comp arison with mic rolensing IGltHGo uld ct al. 2010; SlO- Sumi ct al. 201CJ) , direct imaging (L07-Lafrcnierc ct al. 2007; B15- Bowlc r et al.l 
120151 . and RV (M14 4Montet et al.|l2014l) surveys. The vertical, red lines in each plot indicate the median values for the “hot-start” models 
(solid) and “cold-start” models (dashed). The median values and their 68% confidence intervals are listed in table [l] 


The priors on ujij and are uniform between [0, 2n), 

and that on cosily is uniform between [0,1]. At this 
point, each planet is described by the set of parame¬ 
ters {irip,i,j, aij, These parameters, 

along with properties of the host stars (which we describe 
in the following sections), will be used to determine the 
appropriate observables for each survey. 

4.2. Microlensing Observables 

In order to compare with microlensing surveys, we need 
to map the planet masses and semimajor axes of pop¬ 
ulations of planets into planet-to-star mass ratios and 
projected separations in units of the Einstein radius, 
i.e. , dij) —» ( qi,j,Sij ). This transformation re¬ 

quires that we associate each planet with a host star 
acting as a lens in a microlensing event. We generate 
a ensemble of microlensing events for each planet pop¬ 
ulation by drawing event parameters from mostly naive 
prior distributions. We then calculate posterior distribu¬ 
tions of parameters by weighting each event by its associ¬ 
ated event rate—a statistical probability based on stellar, 
kinematic, and Galactic properties. 

The priors from which we draw event parameters are 
as follows. We choose to fix the source in each microlens¬ 
ing event in the bulge at a distance of = 10 kpc. 

In reality, there is some distribution of source distances, 


however there exists a strong prior on D s . This prior 
is based on two facts: 1) the vast majority of sources 
are in the Galactic bulge, and 2) a majority of lenses 
are also located in the bulge. Since the lens must 
be located in front of the source, the range of D s is 
constrained to be a few kpc, placing most sources behind 
the Galactic center and in front of the back edge of 
the bulge. If we were to allow the source distance to 
vary, it would not significantly change our results, so we 
choose to fix it to simplify our calculations. The lens 
distances, Di^j, are drawn from a uniform distribution 
ranging between [0 kpc,Z? s ] = [0 kpc, 10 kpc]. We draw 
lens masses, from a mass distribution function 

we describe later in this section that includes brown 
dwarf, main sequence, and white dwarf lenses, with 
(initial) masses ranging between [0.01 M 0 ,8 Mg]. 
We choose a random line-of-sight for each event with 
Galatic longitude, h j , and latitude, bi j , drawn uni¬ 
formly between -8°’< l < 0° and -5° < b < 8°, 
respectively. These ranges in l and b are chosen be¬ 
cause they cover the most densely monitored fields 
by the Optical Grav itational Lensing Experiment 
(OGLE; llldalskil l200l . Finally, we determine the 
distribution of lens-source relative proper motions, 
fij j , drawing two components in the Galactic north, 
and east, directions uniformly between 
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Figure 17. Likelihood functions of the parameters of our planet 
population model derived from individual exoplanet surveys across 
the 68% confidence intervals we infer from our order of magnitude 
comparison for each parameter from the combination of all surveys, 
assuming the remaining three parameters are at their inferred me¬ 
dian values. The vertical, grey lines indicate the median values of 
each parameter listed in table [l] for the “Hot-Start” results. Note 
that the o range lines (representi ng the likelihood according to the 
results of ILafreniere et alJI2007h in the upper right and bottom 
panels are included, but are not visible, as they are behind the 
red lines at a likelihood of one, since in each of these panels, a ou t 
is held fixed at its median value, which lies interior to the inner 
working angle of the GDPS. 
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Figure 18. Identical to figure [17] but utilizing the inferred median 
values and 68% confidence intervals inferred for the “Cold-Start” 
comparison results listed in table [l] 


[0 mas yr 1 ,20 mas yr 1 ]. Each planetary microlensing 


event is uniquely described by the set of parameters 

{rifl'p,i,j ■) Oj i^j , i i,j i -Mq, i,j j , Adi i j , , bi^j , } 

The planet-to-star mass ratios for each system are then 
qi,j = m Pt ij/Mijj and we compute the projected sepa¬ 
ration in physical units, r± i ij, by projecting each orbit 
onto the plane of the sk y (see equations 13-20 in § 4 of 
iClanton &; GaudHl2014aD . We then normalize these pro¬ 
jected separations by the Einstein radius, of the 

associated lensing system, Sij = / Re,ij- The Ein¬ 

stein radius is a function of the lens mass and distance, 
as well as the source distance, and is given by 


Re {Mi, Di,D s ) = 


AG Mi D t (D s -DM 1/2 


3.2 AU 


D s 

Mi 


0.5 M 0 
x(l — x ) 11//2 


1/2 


D s 


( 20 ) 

V 1/2 


10 kpc ) 


0.25 


where x = DJD S . 


4.2.1. The Rate Equation 

Microlensing events, are rare, transient phenomena re¬ 
sulting from chance (and very precise) alignments of a 
source star and a foreground lens. The rate of a mi¬ 
crolensing event with given lens mass, Mi, lens distance, 
Di , and lens-source relative proper motion, fi, along a 
given sight line (l,b) is 

d 4 V OD d 2 r dr 
dDidMid^ ~ REV ™ lV ^IdMi ’ (21) 

where Re is the Einstein radius in physical units, u re i 
is the lens-source relative (linear) velocity, v is the local 
number density of lenses at the event location determined 
by Di and (l,b), d 2 T/dfi is the two-dimensional proba¬ 
bility density of a given lens-source relative prop er mo¬ 
tion, fi, a nd dT/dMi is the lens mas s function (see lGriestl 
119911 and lKiraea & Paczvnskilll994l for more information 
on robust calculations of event rates). The lens-source 
relative velocity in physical units is simply v le i = |/x| D/. 

We provide a brief description of these terms to explain 
how relative event rates for an ensemble of microlensing 
events are calculated. For a comp lete explanation, refer 
to § 5 of lClanton fe Gaudil (l2014al) (specifically equations 
33-54). There are a f ew differences in our assum ptions in 
this paper relative to lGlanton fc Gaudl (|2014al lbl). which 
we note below, such as the lens mass function and the 
velocity dispersion of bulge stars we adopt. 


4.2.2. Lens Mass Function 

We adopt a lens mass function consisting of brown 
dwarf, main sequence , and white dw arf lenses similar t o 
those of lGouldl (j200(H and model 1 of lSumi et al.l (120 1 111 . 
We choose not to include contributions from black hole 
and neutron star lenses, as these lenses only constitut e 
a few percent of all microlensing events dGouldl I2000H . 
so their exclusion will not significantly affect our results. 
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Our distribution of “initial lens masses”, Mp, n j t , is thus 


dN 


dlogM^init 


' M 0 - 51 

1V1 l, init 

M~ 0 ' 3 

JW Z,init 

>MrL 


for 0.01 < Mi 5 i n it /Mq < 0.08 
for 0.08 < Mi t i n it /Mq < 0.70 
for 0.70 < Mi ; i n it/ Mq < 8.0. 

( 22 ) 


Our distribution of “final lens masses”, Mi, is identical 
up to Mi < 1.0 Mq, but the white dwarf lenses will have 
final masses described by a Gaussian, dN/dMi | WD oc 
exp [(Mi — AtWD) 2 /( 2 cr| VD )], with a mean mass of pwn = 
0.6 M ( ^ and a st andard deviation of ctwd = 0.16 Mq 
dSumi et al.ll201lT ). Figure [TT)1 plots the initial and final 
lens mass functions we adopt. 

Our final mass distribution yields ratios of the numbers 
of brown dwarf, main sequence, and white dwarf lenses 
of (38:52:10). The fractions of total mass contained in 
each class of objects is (6:69:25), and the fractions of 
microlensing events (weighted by event rate) due to each 
object class is (17:65:17). Thes e frac tions are consistent 
with those reported bv IGouldl (12000D . but note that we 
do not include black hole and n e utron star lenses, and 
the numbers reported by IGouldl • 2000 : group together 
the brown dwarf and main sequence lenses. 


4.2.3. Local Number Density of Lenses 

Since we only need to compute the relative event rates 
for a given ensemble of microlensing events, we can use 
the local mass density, p(x,y,z), rather than number 
density, v(x,y,z), to compute the event rate for any 
particular lens. Here, (x, y, z) are Galactocentric coor¬ 
dinates. We adopt a two-component Galactic model, 
comprised of an axisymmetric (“double-exponential”) 
disk with a 1 kpc hole and a barred, anisotropic bulge. 
T his Galactic mode l is des cribed by equations 39-46 
of iClanton fc Gaudil (l2014af) , including the coordinate 
transformations between (Di,l,b) —» (x,y,z) and the 
normalizations of the bul ge and disk mass density func- 
tions. We demonstrate in lGlanton fc Gaudil f: 2014al ) that 
this Galactic model produces optical depths to microlens¬ 
ing event s that are in rough agreement with those re¬ 
ported by lHan fc Gouldl (gOOl. 

4.2.4. Lens-Source Relative Proper Motion Distribution 

We assume that the velocity dispersion of lenses and 
sources are Gaussian. For lenses in the dislQ, we as¬ 
sume a mean velocity in Galactic north and east direc¬ 
tions of v z = (vp N , vp E ) = (0, u rot ) - (0,10) km s -1 with 
a dispersion of cr v = (o’ v ,n,o' v ,e) = (20,30) km s _1 , 
where u ro t = 220 km s -1 . Lenses/sources that are 
located in the bulge have mean velocities of V;/ s = 
( v i/s,NiVi/ s ^e) = (0,0) km s -1 with a dispersion of er v = 
(<x v ,N, & v ,e) = (130,130) km s -1 . Dividing these veloc¬ 
ities and their dispersions by either Di or D s converts 
them into proper motions. Note that we have adopted a 
slightly larger velocity dispersion for th e objects in the 
bulge than in iGlanton fc Gaudil (|2014al lbl) 130 km s" 1 

versus 100 km s -1 - in order to better reproduce the 
observed d i stribu tion of event timescales measured by 
iSumi et abl (|201lH . We further justify this choice later in 

1 Sources in the disk would have the same mean velocity and 
dispersion, but we assume all sources are located in the bulge. 


this secti on when we descr ibe our comparison with the 
results of ISumi et al.l (1201 If) . 


4.2.5. Statistical Weights for Each Event 

The statistical weight we assign to each microlensing 
event according to the event rate is simply the product 

dr 

I 'Cy = pi,j (x, y, z) — . (23) 

These weights are applied to (probabilistically) correct 
for the fact that many of the stellar, kinematic, and 
Galactic properties are drawn from uninformed priors 
and the fact that we do not know with certainty whether 
or not a given lens is a bulge star or a disk star. These 
weights do not indicate the likelihood of a planet popu¬ 
lation given the results of microlensing (or other) stud¬ 
ies. Note that since we draw our lens masses from an 
informed prior, the events will already be weighted ac¬ 
cording to dr /dMi and we do not need to explicitly in¬ 
clude this factor in F* j. 


4.3. Application to Microlensing Surveys 
4.3.1. I Gould et al\ \201(\ ) 

We apply t he actu al dete ction sensitivities (“triangle 
diagrams”) of lGould et al.l (120101) to our derived distri¬ 
butions of ( qi,j,Sij ) to determine the number of planet 
detections their survey should detect, -/Vd e t,Gio,i and thus 
the planet occurrence rate they would infer, Qi, for each 
population. A likelihood is then assigned to each planet 
population according to t he l ikelihood function we de¬ 
rive for the iGould et alJ (120101) study, .Sfg (plotted in 
figure [3]) , as described in § 13.1.11 

The mathematica l formalism behin d this c alculation, 
developed in § 4 of IGlanton fc Gaudil (l2014al ). is as fol¬ 
lows. We implicitly assume a planet distribution function 
of the form 


d 8 N p \ 


dN p i 

dNpi 

dNpi 

dNpi 

d { 0 } 

fi\ens,i,j 

di 

Lj de 

id duj 

ij dM ° 


d 2 N pi rfiVpi 

d log m p dloga t j dM t 


x 


dN pl 

dDi 




(24) 


such that the joint distribution of log <7 and logs for a 
given realization (i.e. planet population) is 


__ 

dlogq d\ogs dlognip dloga 


Ti 



d {a} 


x TT7 <5(logg-logg iji )<5(logs-logs i j) 

^ 1^1 filens,i,j 

x (5(log m p - \ogm P: ij)5(loga — loga^-) , (25) 


where N ens is the total ensemble number of planets of 
the ith. population. The term J-t/N ens in equation (1251) 
normalizes the joint distribution of log q and log s such 
that there are a total of J'j planets per star in the region 
of log g-log s space we consider. If e ach of the Nq w = 13 
microlensing events of lGould et al.l (|2Q10T ) is indexed by 
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Figure 19. The left panel plots the distribution of “initial lens masses”, Mi i n } t , while the right panel plots the distribution of “final lens 
masses”, Mi. The need to show both distributions arises from the fact that white dwarf lenses will have a mass distribution that is different 
from their initial mass distribution. The relative number of white dwarf lenses is drawn from the initial mass distribution, each of which 
is then given a final mass according to a Gaussian centered at 0.6 Mq with a standard deviation of 0.16 Mq. See text for details. 


then the number of detections we expect for a single 
event in a given realization is then 


N, 


det,G10 ,i,k 


Tv 


d log q J d log s j d log to p J d log a 
d A N p \ 


d\ogq dlogs dlogTOp dloga 
x r ij$(log< 7 ,logs) det ,k 5 


(26) 


where 4>(log q, log s)d e t,fc is the detection limit for the 
fcth microlensin g event, having th e form of the “triangle 
diagrams” of IGould et al.l (I2010D . There is no analytic 
form of <l>(logg,logs)det,fc, so these calculations are per¬ 
formed numerically. Th e to t al nu mber of expected de¬ 
tections for the IGould et al.l (120101) survey for a specific 
planet population is then just the sum over all events, 

iVdet.GiCM = ^Vdet.GiCM.fc, and the expected inferred 
k 

planet occurrence rate is thus Qi = iVd e t,Gio,i/-^Gio- The 
likelih ood of eac h popu lation according to the results of 
the IGould et al.l (120101 ) survey, is given by the 

function plotted in figure [3] We present and discuss these 
results in § [5j 

4.3.2. \Sumi et~al] S201(l) 

We construct the distribution dN p \/ log q\ t for each 
planet population and measu re its sl ope o ver th e range 
of mass ratios to which the iSurni et all (I2010D survey 
is sensitive. Specifically, we bin values of logqy,-, each 
weighted by the associated event rate, Tjj, and use a 
least-squares fit to determine the slope of the resultant 
distribution in the range —4.5 < logg < —2. A likelihood 
is then assigned to each planet popul ation according to 
the likelihood function inferred by the lSumi et al.l (|2010|) 
study, jSf S io (plotted in figure [3j>, as described in 1) 13.1.11 
The results of this comparison are reported in § [5] 


4.4. Direct Imaging Observables 

To determine whether or not a given planet in a 
given population will be detectable by a direct imag¬ 
ing survey, we must derive the mapping () —>• 
(Amag i j, where Amag^ • is the difference in mag¬ 

nitude in a given band (e.g. H or K s ) between the planet 


and its host star and r±,ij is the planet/star projected 
separation in physical units. Since we are testing whether 
or not a single planet population is consistent with all 
surveys, we use the same orbital parameters generated 
for each planet in the comparison with microlensing (see 
§ 14.1.211 to obtain the pr ojection a;., —>■ r \ fol lowing 
equations 13-20 in § 4 of lClanton fc Gaudill2014al 

We require the use of planet evolutionary models to 
determine an absolute (e.g. if—band) magnitude for 
the planet given its mass and age, ■ age pi j). 

This introduces some amount of unquantified uncertainty 
in our analysis, but is a necessary step (see § 16.41 for 
an in-depth description and discussion). We perform 
two separate analyses, one for each class of evolution¬ 
ary models. We use the “cold-start” models at 5x Solar 
(pl anet) atmo spheric metallicity reported in Table 2 of 
iFortnev et al.l (I2008D . To fill in this sparse data, we inter¬ 
polate and extrapolate such that we end up with cooling 
curves with ages log(age/yr) £ [6,10] at uniform sam¬ 
pling of Alog(age/yr) = 0.01, for masses (m p /Mj up ) £ 
[0.5,13] at a uniform sampling of A(m p /Mj up ) = 0.01 
(Fortney, J. J., p rivate commu n icatio n). For the “hot- 
start” models of iBaraffe et al.l ([2003D . we downloaded 
isochrone^] and extrapolated and i nterpolate d in an iden¬ 
tical manner as we did for the IFortnev et all (I2008D mod¬ 
els. 

We then convert the absolute magnitude 
a,ge p i j) derived from the evolutionary 
models to an apparent magnitude, h p ^j(H Pt ij,dij), 
given the distance to the system, D ki j. The difference 
in (apparent) magnitude between the planet and host 
star is then simply 


A hij — H p ij 5 


1 - log ( —hlA 
pc 


— h 


> 


(27) 


where h kt ij is the apparent if—band magnitude of the 
host. 


4.5. Application to Direct Imaging Surveys 
4.5.1. I Bowler et al. I 12015 1) 


2 http://perso.ens-lyon.fr/isabelle.baraffe/C0ND03_models 
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We apply t he ex act same detection limits used by 
iBowler et al.l (120151) in their statistical analysis of the 
PALMS sample to calculate the number of expected de¬ 
tections given a particular planet population (Bowler, B. 
P. private communication). These have the same form as 
the co ntrast curves presented in Table 5 of lBowler et al.l 
(j2015f ) . but are different sets of limits. For each of the 72 
single M stars to which we compare (see § 12.21 for a de¬ 
scription of this sample), we were provided information 
on which filter (either H or K s ) was used for observa¬ 
tions of the target, an apparent H or K s magnitude, 7cr 
contrast limits as a function of angular separation (in 
arcseconds), and FOV coverage fraction as a function of 
angular separation. We obtain distances and ages (com¬ 
plete with un certai nties) for each star from Table 2 of 
IBowler et al.l (120151) . 

For a given planet population, *, we determine the 
number of expected detections around each of the 72 
single M stars in the PALMS sample, which we will enu¬ 
merate by the subscript k, in the following manner. For 
each planet, j, in this particular popula tion we draw a n 
age from the estimates provided bv IBowler et al.l ( 2015 ). 
For stars in young moving groups, we assume Gaussian 
uncertainties, while those with an age range listed, we 
assume a distribution that is uniform in log-space. Us¬ 
ing the age for each planet around a given star we draw 
from one of these distributions, age i j>fe , the distance to 
that system, D*^, and the apparent magnitude of the 
star in the appropriate band, we determine the planet-to- 
star contrast Amag^ k as described above. Combining 
this contrast with the projected separation we calculate 
for each system, r±^j, the FOV coverage fraction at this 
projected separation, and the stability limits for the host, 
Ostab,fc, (if applicable, see below for description), we de¬ 
termine whether or not such a planet is expected to be 
detected by the PALMS survey. 

We developed a generalize d mathematical formalism 
for such calculations in § 4 of lClanton fe Gaudil (|2014al) . 
which we apply here. We implicitly assume a planet dis¬ 
tribution function of the form 


d 9 N p i 
d {a} 


dNpi 

dNp\ 

dN p \ 

dNp\ 

di 

i,j de 

id du} 

ij dM ° 


x 


d 2 N p i 

dNp\ 

dV pl 

dlogm p dloga 

i,j d{ mag*) 

id dD * 


dN p i 

d{ age) itj 


(28) 


where mag* is the apparent magnitude of the hosts in 
either H— or K s — band, such that the joint distribution 
of Amag and r_L for a given realization (i.e. planet pop¬ 
ulation) is 


d 4 W 


pi 


d(Arnag) dr± d\ogm p dloga 
d 9 N p i 


Ti 


N t 


ens J {a} 


d {a} 


d {a} 


5(Amag - Amag, )d(r ± - r± ti<j ) 


x 6(logm p - log m Pt i t j)6(log a — log a itj ) , (29) 


normalizes the joint distribution of Amag and r± such 
that there are a total of J-) planets per star in the region 
of Amag-ri space we consider. The number of detec¬ 
tions we expect for a single target in the PALMS sample 
(and in a given realization) is then 


Adet,Bi 5 ,i,fc = J d( Amag) J dr± J d log m p J dlog< 

d^Npi 


d(Amag) dr± dlogm p dloga . 

X $(Amag,rj_)det,fc^(r_L)FOV,fc^(a) S tab,i: , 

(30) 


where <E>(Amag, rj_)det,fc is the contrast limit, <h(rj_)Fov,fc 
is the FOV coverage fraction as a function of r± (a num¬ 
ber between zero and one), and < f>(a) s t a b,fc = 0(cistab,fe — 
a) is the planet stability criterion as a function of semi¬ 
major axis for the fcth star in the PALMS sample pro¬ 
vided to us (Bowler, B. P. private communication). 
Some of the stars in the PALMS sample have known 
(wide-separation) binary companions that are expected 
to threaten the long-term stability of some p l anetar y or¬ 
bits about the primary (Holman&Wlegei't I1999D . In 
their statistical analysis. IBowler et alJ ( 20151) assume no 
sensitivity to planets with separations that would place 
them in unstable orbits according to t he st ability crite¬ 
rion developed in iHolman fe Wiegertl (I1999D . We adopt 
these same stability limits and assume there is no sen¬ 
sitivity to planets on orbits beyond the maximum sta¬ 
b le se mimajor axis, a s t a b fe, calculated bv IBowler et all 
(120151) . 


There are no analytic forms of 4>(log(7,logs)det,fe and 
< h( r _L)FOV,fc, so these calculations are performed numer¬ 
ically. The total number of expected detections for the 
PALMS survey given a specific planet population is then 


just the sum over all targets, Vdet.Bis.i = Vdet.Bis.i.fc- 

k 

The likelihood of each population according to the re¬ 
sults of the PALMS survey, which reported zero plane¬ 
tary companion detections over all their targets, is then 
det,Bi 5 ,i) = exp(- V det ,Bi5,i)- This process is re¬ 
peated for many different planet populations to build up 
constraints on the parameters of our assumed planet dis¬ 
tribution function, which we present and discuss in § [5] 


4.5.2. I Lafreniere et al\ K200 

We determine the number of expected detections for 
the GDPS given a particular planet population in a sim¬ 
ilar fashion as we describe above for the PALMS sur¬ 
vey, with a few minor differences. For each of the 16 M 
stars in the GDPS, we obtain H —band mag n itudes, dis¬ 
tances , and age ranges from Table 1 of lLafreniere et all 
• 2007 :, and we obtain 5 a detection limi ts from thei r Ta¬ 
ble 4. The detection limits provided by [Lafrenier e et all 
(120071 ) are magnitude differences in the NIRI CH4-short 
filter (1.54-1.65p,m) as a function of projected separa¬ 
tion (in arcseconds). For planets in a particular popu¬ 
lation, we compute the expected contrast in this CH4 
filter, AracH 4 ,ij , as 


where N ens is the total ensemble number of planets of 
the ith population. The term fFi/N ens in equation (1291) 


Am C H4,i,i = - 2-5 log > (31) 
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where Ap ? is the magnitude difference in if—band given 
by equation (EH, /ch 4 and fjj are the mean flux den¬ 
sities of the planet in the NIRI CH4-short and H fil¬ 
ters, respectively, and the host star magnitudes in these 
two filters are assumed to be equal (see equation 2 
and surrounding text in § 4.1 of lLafreniere et ahll2007lf . 
Since the ratio of flux densities in these filters typically 
lies in the range 1.5-2.5 for giant pl anets, we adopt a 
value of fnru / f h = 2 for all planets (lAllard et all 120011 : 
IBaraffe et aiJl2005 lLafreniere et al.ll2007f) . 

The mathematical formalism behind this calculation is 
identical to that presented for the PALMS survey above, 
but with $(rj_) F ov,fc -» 0Jr_L - (pFOV,L07-D*,i,j)] and 
c f > (a)stab,fc -> 1 in equation ([30]), where pfov,lo7 = 10" 
is the (radial) FOV of the GDPS, out to which we as¬ 
sume they have a unifo r m cov erage fraction of one. We 
follow lLafreniere et al.l (1200711 and do not adopt a sta- 
bility criterion for bin aries in our comparison (note that 
lLafreniere et al.ll2007l only report one of their M stars as 
having a known companion, GJ 234AB). 

Thus, for each planet population, we numerically com¬ 
pute a number of expected detections for the GDPS, 
Adet,L07,i,fc) where the subscript k enumerates the 16 
M stars in the sample. The total expected detec¬ 
tions for their survey (for a given population) is then 

Adet,L 07 ,i = ^fV d et,L 07 ,i,fc- We assign a likelihood 

k 

to each population according to the likelihood function 
-%07(^Vdet,L07,i) = exp(-7V d et,L07,i) since the GDPS dis¬ 
covered no planetary companions to any of their targets. 
Repeating this process for many different planet popula¬ 
tions, we derive the constraints presented and discussed 
in §03 

4.6. RV-Trend Observables and Application to the CPS 
TRENDS Survey 

We compute the number of expected long-term RV 
trend detections for the CPS TRENDS survey given 
a specific planet population by calculating the magni¬ 
tude of acceleration of the host star due to a planet 
on a given orbit and comparing with the detection sen¬ 
sitivities of the sample. Again, since we are testing 
whether or not a single planet population is consis¬ 
tent with the results of all the surveys we consider, we 
use the same orbital parameters we generated for the 
comparison with microlensing and direct imaging (see 
§ 14. 1.21) . This allows us to immediately compute the ex¬ 
pected acceleration for each planet in a given population, 

using equation m- We also 
need to know the periods of the planets’ orbits, Pij, to 
determine their detectability as a long-term trend, rather 
than as “secure” detections. For a planet to be detected 
as a long-term RV trend, the period of the planet must 
be larger than the baseline of observations, such that the 
RVs of the host ar e monitored for only a fr action of an 
orbital period fsee lClanton fc Gaudill2014al lb[). The or¬ 
bital period of each planet is easily computed via Kepler’s 
Third Law (equation [TT]) given the masses of the CPS 
TRENDS M dwarfs and the planet’s semimajor axis. 

4.6.1. I Montet et al . I H201A) 

For each star in the CPS TRENDS survey, which 
we will enumerate with the index k, we determine 


the “threshold” acceleration as a function of orbital 
period, above which planets are expected to be de¬ 
tected as long-term RV trends, 7>thr ,k(P), directly from 
equation m■ We describe the reasoning behind 
this minimum-detectable acceleration in § 13.1.31 and in 
iClanton fc Gaudil (j2014al [bh . Previously, in the order of 
magnitude comparison presented in § 13.1.31 we used the 
median sample parameters of the CPS TRENDS survey 
to derive rough constraints on the properties of planets 
consistent with the survey results. Now, we will per¬ 
form a more accurate calculation, including a different 
detection limit for each star in the sample that is de¬ 
termined from the specifics of the observations of each 
target. We obtain masses, M*^, for each of the 111 M 
dwarfs monitored by the CP S TRENDS survey from Ta¬ 
ble 1 of lMontet et al.l (I2014f) . and from their Table 2, we 
obtain the total number of observations, Nk, the base¬ 
line of observations, T*,, and the total RV measurement 
uncertainty (the quadrature sum of the estimated pho¬ 
ton noise and stellar jitter), ak, for each target. Thus, 
for each star, we compute the threshold acceleration (i.e. 
trend detection limit) as a function of orbital period us¬ 
ing equation Cl to obtain the “phase-averaged” S/N, 
Qk, and substituting this into equation (1151) . An exam¬ 
ple of such a detection limit is shown in figure [9] 

We implicitly assume a planet distribution function of 
the form 


d 7 N p i 
d {a} 


R y,i,j 


dN pl 

dN p i 

dNpi 

dNp\ 

di 

id de 

id du 

. . dMn 

l,J u 


d 2 N , 


pi 


dlogm p dloga 


dN r 


pi 




dM* 


(32) 


such that the joint distribution of v and log P for a given 
realization (i.e. planet population) is 
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where N ens is the total ensemble number of planets of 
the ith planet population. The term Pi/N ens in equa¬ 
tion (1331) normalizes the joint distribution of v and log P 
such that there are a total of T r planets per star in the 
region of v-log P space we consider. The number of long¬ 
term RV trend detections we expect for a single target in 
the CPS TRENDS sample for a given population is then 
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(34) 


where $(i>) tr ,fc = 0[fi(logP) -'Cthr.fe(logP)] is the ac¬ 
celeration detection limit for the fcth target at the or¬ 
bital period of a given planet and where <h(log P)t r ,fc = 
0(logP — logTfc). The first detection criterion simply 
requires that the acceleration of the host star due to 
a planet has large enough magnitude to be detectable, 
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while the second criterion ensures that the period of the 
planetary orbit is longer than the baseline of observa¬ 
tions. The total number of expected long-term RV trend 
detections for the CPS TRENDS survey for a specific 
planet population is then just the sum over all targets 
in the sample, jV tr .,; = and the likelihood of 

k 

each planet population is given by the likelihood func¬ 
tion Jz= (iVtr,i) 4 exp(—-/V trj j)/4!, which we plot 
in figure [3] We present and discuss the results of our 
comparison with the CPS TRENDS survey in § [5] 

5. RESULTS AND DISCUSSION 

We find that the results of microlensing, RV, and di¬ 
rect imaging surveys are consistent with a single planet 
population described by a joint power-law distribution 
function in planet mass and semimajor axis. We present 
our final results in figures [23 (“hot-start”), UHl (“cold- 
start”) and m and summarize the median values and 
68% confidence intervals for each of the four parame¬ 
ters of our planetary distribution function, which are 
consistent with all surveys we consider, in tabic [5] We 
plot the associated likelihoods of the individual surveys 
across these 68% confidence intervals in figures [28] and 
1291 demonstrating that the resultant population of plan¬ 
ets has a reasonably high likelihood of producing the ob¬ 
servations reported by each survey. 

There are several subtleties embedded in the con¬ 
straints obtained from the synthesis of all survey results, 
so in order to demonstrate that we understand the re¬ 
sults of our calculations, we first present the constraints 
obtained from each individual survey and compare with 
those we found in the order of magnitude evaluations of 
§ [3] before reporting our final results. 

5.1. Constraints from Individual Surveys 
5.1.1. Microlensing 

Gould et al. (|2010f ): Figure [20] shows the likelihood 
contours on the parameters of our planet distribution 
function that are consistent with the results of this sur¬ 
vey. These constraints have some differences with those 
we obtained in the order of magnitude evaluations, but 
the general intuition we develop in S I3.1.1l for understand¬ 
ing these constraints still holds. 

As we found in § 13.1.11 the most important constraint 
here is on the combination of a and A shown in panel (b) 
for the same reasons. We find that the constraints on the 
combination of /3 and A shown in panel (d) are more pro¬ 
nounced because, by applying different det ection limits 
for ea ch of the 13 microlensing events in thc lGould et alJ 
(2010) sample, we are effectively sampling the planet dis¬ 
tribution function over various ranges of s (and thus a), 
thus providing more sensitivity to /?. As a result, the con¬ 
straint on the combination of a and ft is slightly stronger 
in this calculation relative to those of § 13.1.11 

Another difference here is that there is no minimum 
value of ao U t (within the range we consider) as we found 
for the order of magnitude calculations. There, we had a 
sharp cutoff at 1.4 AU due to the sensitivity criteria we 
applied, which is not apparent here. However, panel (a) 
shows that small values of a ou t req uire very large nor - 
malizations to be consistent with the lGould et al.l ([2010( 1 
results. Finally, in panel (c), we see that very large and 



Figure 20. Likelihood contours as a function of pairs of pa¬ 
rameters describing planet distribution functions that are found 
to be consistent wi th the measuremen t of the integrated planet 
frequency, Q, by the lGould et al.l (|2010T ) microlensing survey. Con¬ 
tours are drawn at levels of 68% and 95% of the peak likelihood 
and are marginalized over all other parameters except those being 
plotted. 

v ery small valu e s of a are excluded unless a ou t is small. 

Sumi et al. (120101) : Figure [21] shows the likelihood 
contours on the parameters of our planet distribution 
function that are consistent with the measurement of 
the slope of the mas s ratio dist r ibutio n function (p = 
—0.68± 0.20) by the iSumi et al.l (12010H study. It is ob¬ 
vious why the only constraining power offered by this 
survey is on a, but we note that the contours are not 
centered at a = —0.68 and the 68% contours spread 
larger than 0.20. Since q = m p /Mi, the conflation of 
dNpi/d\ogm p with the lens mass function (see 14. 2. 2D 
means that dN p \/d\ogm p dN p \/d log q. This is in con¬ 
trast to the order of magnitude constraints presented in 
§ra where we had assumed a characteristic lens (host) 
mass rather than marginalizing over the population. 

5.1.2. Direct Imaging 

iBowler et al.l (1 201 .51 ) : We plot the likelihood con¬ 
tours derived from our comparison with the results of 
the PALMS survey in figure [22] (see § 14.5.11 for method¬ 
ology). These constraints are very similar, qualitatively, 
with those derived in our order of magnitude analysis of 

§nn 

The PALMS survey can only weakly constrain combi¬ 
nations of a ou t with a, /?, and A. As shown in panels (c) 
and (e), large values of j3 and a are disallowed for large 
values of a ou t, since we would expect too many direct 
detections from such populations. In panel (a), we see 
the same behavior with large values of a ou t and A, al¬ 
though here it is less pronounced since this set of contours 
marginalizes over all values of a and fd. For sufficiently 
small values of a and/or /3, most combinations of a out 
and A are allowed. Our results are more sensitive to (d, 
and especially more sensitive to a, than they are to A. 

The major difference between the order of magnitude 
comparison presented in 03.1.2l and the results of our de- 
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Figure 21. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent with the measu rement of th e slope of the mass-ratio dis¬ 
tribution function by the ISurni et al.l l|201 [j') microlensing survey. 
Contours are drawn at levels of 68% and 95% of the peak likeli¬ 
hood and are marginalized over all other parameters except those 
being plotted. 

tailed calculations we report here is this enhanced sen¬ 
sitivity to a shown in panels (b), (c), and (f). In our 
order of magnitude comparison, we assumed that each 
star in the PALMS sample had the same detection limit 
(the median of the sample), whereas here, we employ 
the individual sensitivities for each star. This allows 
us to sample the planet distribution function at differ¬ 
ent contrasts (and thus, through the use of evolutionary 
models, different planet masses), providing more sensi¬ 
tivity to the value of a. Not surprisingly, this is also 
where we see the largest difference between the “hot-” 
and “cold-start” models (see panel c). However, as we 
will demonstrate in § 15.21 there is not a significant dif¬ 
ference in the constraints derived from the “hot-” and 
“cold-start” models over the regions of parameter space 
al lowed by a .joint compariso n with all surveys. 

lLafreniere et aD (|2007f) : Figure 1251 plots the like¬ 
lihood contours on the parameters of populations con¬ 
sistent with the GDPS we derive following the proce¬ 
dure outlined in § 14.5.21 The interpretation of these 
constraints is identical to that discussed above for the 
PALMS survey. Although the GDPS is, in general, a 
deeper survey than PALMS, the large difference in sam¬ 
ple sizes between the surveys means that the GDPS con¬ 
straints are weaker than those derived from the PALMS 
survey. 

5.1.3. Long-Term RV Trends 

Montet et al. (l2014h : We plot the likelihood con¬ 
tours on the parameters of our planet population model 
that are consistent with the CPS TRENDS survey in fig¬ 
ure [U (methodology for this comparison is presented in 
§ 14.6.11) . The constraints are qualitatively very similar 
to those of our order of magnitude comparison shown in 
§ 13.1.31 As expected, the detection of RV trends places 
a lower bound on a ou t shown in panels (a), (c), and (e). 


Figure 22. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent with the non-detection o f any planetary com panions in 
the PALMS direct imaging survey <\Bowler et al.|[2015lh Contours 
are drawn at levels of 68% and 95% of the peak likelihood and are 
marginalized over all other parameters except those being plotted. 
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Figure 23. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to be 
consistent with the non-dete ctio n of any planetary companions in 
the GDPS HLafreniere et afll2007l) . Contours are drawn at levels 
of 68% and 95% of the peak likelihood and are marginalized over 
all other parameters except those being plotted. 

We find an anticorrelation between a and A (panel b) 
such that smaller values of a, for which more planets are 
distributed to lower masses and thus induce smaller ac¬ 
celerations on their hosts, require larger normalizations 
in order to be consistent with the observed number of 
trends. The CPS TRENDS survey places a constraint 
on the combination of a and (3 (panel f); larger values 
of /3 distribute more planets to larger separations where 























































































































































25 


11 
O) 

^ 0 


\ 

(a) 

'— 



Constraints: 

68% Contours Montet+ < 2014 > 

95% Contours 



Figure 24. Likelihood contours as a function of pairs of param¬ 
eters describing planet distribution functions that are found to 
be consistent with the detection of four significant, lon g-term RV 
trends in the CPS M dwarf sample IMontet et al.1120 141) . Contours 
are drawn at levels of 68% and 95% of the peak likelihood and are 
marginalized over all other parameters except those being plotted. 

they might be detectable as trends, and thus smaller val¬ 
ues of a are required to keep the number of planets mas¬ 
sive enough to induce detectable accelerations consistent 
with the observations. 

5.2. Synthesized Constraints 

We have shown that, individually, the results of any 
single exoplanet discovery survey can only constrain com¬ 
binations of the parameters of our assumed planet distri¬ 
bution function. However, a joint comparison with the 
results of multiple detection methods enables us to de¬ 
rive meaningful constraints on the individual parameters 
of a planet population that is consistent with the results 
of all surveys we consider. We plot the likelihood con¬ 
tours on pairs of parameters describing a population that 
is consistent with the results of various combinations of 
surveys in figure [25] for our “hot-start” calculations and 


figure [26 


Since the 


for our “cold-start” calculations. 

survey con s trains a com- 


Gould et al.l 


bination of a and A and the ISumi et al.1 (|2010f) survey 
constrains a, the combination of these two surveys pro¬ 
vides a relatively tight constraint on A. With a and A 
more or less fixed, the CPS TRENDS survey places lower 
bound on fd and a ou t, below which there are not enough 
planets at large enough separations to be detected as 
long-term drifts. However, due to degeneracies between 
fd and aout, we are unable to place an upper bound on 
either of these parameters. Very large values of fd are al¬ 
lowed by the trends as long as a ou t is relatively small. In¬ 
clusion of the direct imaging results does not improve the 
constraints much for similar reasons. Arbitrarily large 
values of fd are allowed by direct imaging surveys as long 
as a ou t < 10 AU. This is due to the fact that there is 
a sharp cutoff in detectability at projected separations 
approaching the inner working angle. 

There is not much difference between the results of the 


“hot-” and “cold-start” calculations. The “cold-start” 
calculations are slightly less constraining, particularly on 
aout) but the difference is not much and well within the 
uncertainties (see table [2]). One might assume that this 
is because the difference in expected planet luminosity 
between the two classes of evolutionary models is great¬ 
est at young ages, a region of parameter space occupied 
by just a fraction of the targets in direct imaging sur¬ 
veys. However, this is not actually the case, because it 
is these very targets (i.e. stars with ages < 100 Myr) 
that provide the best constraints from imaging surveys 
as overall; most planets are not expected to be detected 
around the older stars. 

Instead, the lack of a significant difference between the 
“hot-” and “cold-start” predictions can be attributed to: 
1) the difference in the predicted luminosities between 
the two classes of models decreases with planet mass at 
all ages, 2) our assumed power-law form for the planet 
distribution function, and 3) the most likely range of a ou t 
we find lies within a region where direct imaging surveys 
are mostly in sensitive to p la nets . The bottom-left panel 
in Figure 1 of lFortnev et al.l (1200811 demonstrates that, at 
all ages, the predicted luminosities for lower-mass plan¬ 
ets differ less between the “hot-” and “cold-start” mod¬ 
els relative to the more massive planets. This is impor¬ 
tant, since we find that a is relatively well-constrained 
at negative values, and thus the mass function is heavily 
weighted towards lower-mass planets. Additionally, since 
our assumed planet distribution function has a power-law 
form, small changes in a, fd, and/or a ou t can have a dra¬ 
matic effect on the number of expected planet detections. 
As for the third reason, the detectability of exoplanets 
by direct imaging surveys dramatically decreases as the 
projected separatio n approaches the inner working angle 
(see e.g. figure 2 o f lLafreniere et al1l2007l and figure 3 of 
iBowler et al.ll2015h , which means that we do not predict 
any exoplanet detections with separations a < 10 AU for 
either survey, despite our choice of evolutionary model 
(see figures [721 and l23l). 

We plot the corresponding one-dimensional likelihood 
functions on each parameter, marginalized over the other 
three parameters, in figure 1271 and summarize the me¬ 
dian values (and their 68% confidence intervals) in ta¬ 
ble [2] for both “hot-” and “cold-start” planet evolution¬ 
ary models. The final parameter constraints we report 
in table [2] do differ from those we found in our order of 
magnitude comparisons (presented in table[l]). While the 
median values are consistent between the two indepen¬ 
dent analyses, the results of our detailed calculations are 
more constraining, particularly on fd and a ou t (compare 
figures [16] and [27]). 

The likelihood distribution of a we find in our detailed 
analysis is slightly shifted towards smaller values, but 
has a 68% confidence interval with no larger or smaller a 
“spread” than that derived from our order of magnitude 
analysis. The same is true of A, although the shift is 
less pronounced. Although consistent between the two 
analyses, the resultant median value of fd of our detailed 
calculations is a factor ss 2 larger than that of the or¬ 
der of magnitude comparison. Increased sensitivity to 
fd due to the inclusion of (real) variations in detection 
limits within the individual surveys, which thus sample 
our planet distribution function at different semimajor 
axes, means we also find that fd is better constrained 
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Figure 25. Likelihood contours as a functi on of pairs of para meter s describing plan et distribution functio ns that are found to b e con¬ 
sistent with the results of microle nsing (G10 -Gould ct al. 20101 : SlO dSumi et all 120100 . direct imaging fL07- fLafreniere et al.l [20071 : B15— 
IBowler et al.1120151). and RV (M 14 4Montet et al.l[2014T ) surveys. Here, we have utilized the predictions of the “hot-start” planet evolutionary 
models of IBaraffce et all (120030 . Contours are drawn at levels of 68% and 95% of the peak likelihood and are marginalized over all other 
parameters except those being plotted. 


relative to the results of our order of magnitude calcula¬ 
tions (see §ED for more discussion). In turn, due to the 
previously discussed degeneracy between /3 and a ou t (as 
well as variations in specific detection sensitivities), we 
find that a ou t is also better constrained by our detailed 
comparison, although the median values between the two 
analyses are indistinguishable. 

We have demonstrated that there is overlap among the 
constraints we derive from a comparison with each of the 
individual surveys, which means that the results of all the 
surveys we consider can be explained by a single planet 
population described by a simple, joint power-law distri¬ 
bution function given by equation ©• As before in our 
order of magnitude comparison, we now determine how 
the likelihood of such a planet population, according to 
each survey individually, varies over the 68% confidence 
intervals we derive on a , /3, A , and a ou t- 


In figures E51 (“hot-start”) and HZ 7 )! (“cold-start”), we 
plot the normalized likelihood functions (i.e. the likeli¬ 
hoods normalized by the maximum likelihood value for 
each individual survey) of given parameters across their 
68% confidence intervals, assuming the remaining three 
parameters are at their median values (listed in table©. 
These figures show that the likelihoods of these param¬ 
eters from each individual survey are consistent to well 
within ler (i.e. their normalized likelihood is greater than 
0.32) across much of their inferred 68% confidence in¬ 
tervals, and have normalized likelihoods > 0.7 at the me¬ 
dian value for all surveys (indicated by the vertical, grey 
lines). This demonstrates that this overlap in parame¬ 
ter space we identify is not improbable, and has a rather 
good chance of being real. Note that while in each of the 
panels there are some regions where the relative likeli¬ 
hood for a given survey seems low (e.g. in the upper left 






















































27 



68% Contours 
95% Contours 

Cold Start 


2 [ 

8 0 

-2 


(b) 


i ... i ... i ... i 



(C) 





2 

SM-W 

ca 0 




-2 

. / / 

// ^ / 

' ' ' * 

* / / 


✓ / 

/ / 


- 2-10 1 2 

log(A [dex 2 ]) 


Constraints: 

--- G10 + S10 
+ L07 + B15 

. G10 + S10 

+ M14 

- G10 + S10 

+ L07 + B15 
+ M14 



0 1 2 
log(a out [AU]) 


-2 


0 

a 


Figure 26. Identical to figure [25] but utilizing the predictions of the “cold-start” planet evolutionary models of IFortnev et all 120081! . 


Table 2 

Median values and 68% uncer tainties on the para mete rs of our planet di stribution mo dels for our detailed c omparisons wi th the results of 
the microlensing surveys by JGould et all <201(11 and lSumi et all <201(11 . the GDPS lll.a.f'reni ere et al.||2007|l and PALMS < Bowler et all 
1201511 direct imaging surveys, and the CPS TRENDS survey llMontet et al.1120141 . 


Planet Evolutionary 
Models 

Median 

Q 

Values anc 

p 

i 68% Uncer 
A [dex -2 ] 

ainties 
&out [AU] 

“Hot-Start” 
(Baraffe et al. 2003) 

—0.86±g;?| 

1 1+ 1 - 9 

0 21 +°- 20 
U - Z1 -0.15 

Iq+26 

iU -4.7 

^^Tlold- Start"" 

(Fortnev et al. 2008) 

-0-85tg;S 

1 1+ 1 * 9 

0 21 + 0 - 20 
U - Zi -0.15 

-,9+50 
±z —6.2 


panel of figure l28l the l ikelihoo d of q 01 it. according to the 
PALMS survey of lBowler et al .1120 151 drops below 0.32 be¬ 
yond log clout ~ 1-2), this does not necessarily imply that 
a portion of these 68% confidence intervals are relatively 
unlikely- keep in mind that the other three parameters 
are being held fixed at their median values in each panel. 
Although these likelihood curves are a little noisy, they 
are remarkably similar (qualitatively) to those we found 
in our order of magnitude analysis shown in figures 1171 


andfTBl 

6. UNCERTAINTIES 

The main sources of quantified uncertainty in this 
study are those of the measurements made by the six 
different surveys with which we jointly compare. By 
employing the full likelihood functions (or approxima¬ 
tions thereof) of these measurements, rather than the 
maximum likelihood values, the uncertainties we derive 
on our results naturally include the contributions from 
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Figure 27. Relative li kelihoods for each para meter, marginaliz ed over the other thre e parameters, derived from our detailed compa rison 
with micr olensing (GlCHGould et al.ll20KJ : SlC HSumi et al.ll20 1 (l i . direct imaging (L07- TCafreniere et al.|[2007l : B15- lBowler et aI1l2015ll . and 
RV I'M 14 iMontet et" Ti l . 2014) surveys. The vertical, red lines in each plot indicate the median values for the “hot-start” models (solid) and 
“cold-start” models (dashed). The median values and their 68% confidence intervals are listed in table [2] 


each individual measurement. However, there are mul¬ 
tiple sources of unquantifiablc uncertainty in this study, 
in particular in the adopted models and priors necessary 
to derive the mapping between the properties of planets 
( m p and a) into the observables relevant to the differ¬ 
ent types of surveys. There are certainly differences in 
the properties of the samples of the individual surveys 
(e.g. distributions in stellar mass and metallicity), the 
magnitudes of which are unknown due to a lack of the 
appropriate observational data. We discuss the nature of 
these unquantified uncertainties below, arguing that at 
least some of them should be subdominant. 

6.1. Priors on Unknown Stellar, Galactic, and 
Planetary Orbital Properties 

As we have demonstrated, a robust, joint comparison 
between the results of these different techniques is con¬ 
tingent on an ability to map basic planet properties into 
the observables relevant to each method. This requires 
knowledge of the physical properties of the hosts as well 
as the planetary orbital parameters, since each technique 
is sensitive to a unique combination thereof. While we at¬ 
tempt to account for uncertainties in these properties by 
marginalizing over prior distributions, imperfect knowl¬ 
edge of their actual distributions leaves some level of 
residual error. 

The routine observables in a typical planetary mi¬ 


crolensing event contain little information about the 
physical properties of the lens system beyond the planet- 
to-star mass ratio, q , and the projected separation in 
units of the Einstein radius, s. Measurement of the pri¬ 
mary mass, and thus m p , requires detection of higher 
order effects in the light curve (e.g. finite source effects, 
microlens parallax, xallarap) to break the degeneracy be¬ 
tween lens mass and distance. The mapping necessary 
to compare our planet population models to the results 
of microlensing surveys, ( m p ,a ) —>• (q,s), thus depends 
on a Galactic model, priors on the lens mass and dis¬ 
tance distribution functions, and priors on the planetary 
orbital parameters. 

In this paper, we adopt the same Galactic model and 
prio r on lens distances as in iClanton fe Gaudil (|2014afl 
and IClanton fe Gaudil (|2014blh We demonstrated in 
those studies that we are able to reproduce the observed 
di stributions o f host star parameters (Re, tE, \u\, Af;) 
of lGould et al.l (120101 ) microlensing survey when weight¬ 
ing them by their expected event rates, with the possible 
exception of the distribution of lens distances. We found 
that events in our simulated sample typic ally have larger 
values of Di relative to the IGould ct al.l (I2010D sample. 
However, given the small number of actual distance mea¬ 
surements in the microlensing samples, this ultimately 
may not be significant. Although we adopt a different 




















































29 





Figure 28. Likelihood functions of the parameters of our planet 
population model derived from individual exoplanet surveys across 
the 68% confidence intervals we infer for each parameter from the 
combination of all surveys, assuming the remaining three param¬ 
eters are at their inferred median values. The vertical, grey lines 
indicate the median values of each parameter listed in table [2] for 
the “Hot-Start” results. 
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Figure 29. Identical to figure l28l but utilizing the inferred median 
values and 68% confidence intervals inferred for the “Cold-Start” 
comparison results listed in table l2l 

lens ma ss function in this study than in lClanton fe Gaudil 
we st ill find agreement b etween our simulated 
events and the iGould et all (12010H sample. We there¬ 
fore expect uncertainties due to our choice of priors and 
Galactic model to be subdominant. 

The priors over which we marginalize the planetary 
orbital parameters (e, i, w, Q, M$) are those predicted 
for randomly-oriented orbits at random phases. We are 
effectively assuming there is no preference of orbital ori¬ 
entation to location in the Galaxy or host star prop¬ 


erties, which would mean that any uncertainties intro¬ 
duced into our analysis by unknown orbital parameters 
should be naturally included in our analysis by virtue of 
marginalizing over these priors. However, we note that 
orbital eccentricity might actually depend on location in 
the Galaxy and/or on host star properties, resulting in a 
potential source of unquantified uncertainty in our anal¬ 
ysis. 


6.2. Stellar Metallicities 

The metallicities of the microlensing and direct imag¬ 
ing samples are unknown, however it is reasonable to 
expect that at least the microlensing sample should have 
a significantly different distribution than those of the RV 
and direct imaging surveys. The RV and direct imaging 
samples are limited to nearby M stars (within tens of 
parsecs), whereas microlensing probes stellar hosts sev¬ 
eral kpc into the Galaxy. Multiple studi es have mea¬ 
sured a Galactic metallic ity gradient (e.g. iCheng et all 
I2012t IHavden et a.1.1120 Mil with metallicity decreasing as 
distance from the Galactic center increases. The metal¬ 
licities of the disk lenses (and at l east a fraction of 
bulge lenses; see iBensbv et al.11201311 are therefore ex¬ 
pected to be enhanced relative to the local RV and di¬ 
rect imaging samples. Despite this, and the fact that 
there is a strong correlation between giant planet fre¬ 
quency and host star metallicity measured by RV surveys 
(lGonzalezlll997t IF scher fc Valentil 120051: I.Tohnson etak 

l2010HMontet et al.H2014li 7 we found in lClanton fc Gaiidl 

( 2014al) that the measurements of giant planet frequency 
made by RV and microlensing surveys are consistent 
without a need to correct for any meta lli city effects. 
This is discussed at length in § 6.4 of lClanton fc GauTIil 
(1201 4all . 

6.3. Non-M Dwarf Lens Contamination in Microlensing 

Samples 

The rate of microlensing events along a given sight line 
depends explicitly on the lens mass distribution function 

since T oc 0e oc (see e.g. lGaudill2012t l. The slope 

of the mass function is such that there are greater num¬ 
ber of low-mass stars relative to higher-mass stars, and 
therefore a majority of lenses are M dwarfs. However, 
there is some level of contamination of the microlens¬ 
ing samples by earlier spectral types (mostly GK), and 
stellar remnants. An estimated ~ 20% of microlensing 
events are due to remnants (primarily white dwarfs) and 
are essentially completely indistinguishable fro m other 
types of lenses by their timescale distribution i Gouldl 
1200011 . Nevertheless, we expect the uncertainties due to 
any such contamination to be small in comparison with 
the Poisson errors due to the small number of planet de¬ 
tections in the microlensing samples that are included in 
our analysis. 

6.4. Planet Evolutionary Models 

For the purposes of comparing our planet popula¬ 
tion models with the results of direct imaging surveys, 
we must derive the expected planet-to-star constrast, 
Amag(m p , age, distance), for each system. Expected ab¬ 
solute H and/or K s band magnitudes of each planet 
are estimated using the planet evolutionary models of 
iBaraffe et al.1 (|2003h and iFortnev et al.l ([2008 1 given the 
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planet’s mass and age; the latter is assumed to be the 
same as the host’s age. These are then converted into 
apparent magnitudes using the distance of the host. The 
contrast is then simply the difference between the appar¬ 
ent magnitude of the planet and that of the host in the 
appropriate band. While the uncertainties in the ages 
and distances of these systems are somewhat taken into 
account in our analysis by marginalizing over the ap¬ 
propriate probability density functions of these measure¬ 
ments (see below for additional discussion), uncertainties 
in the planet evolutionary models themselves are not. 

There are two clas ses of these models, “c old start” (e.g. 
Fortnev et al.l l2008h and “hot start” (e.g. iBaraffe et all 
200311 . which are broadly distinguished by their assump¬ 
tions about the initial conditions of a newly-formed pro¬ 
toplanet concomitant with differe nt formation mecha¬ 
nisms. The “cold-start” models of lFortnev et al.1 (|2008H 
we adopt in this paper follow the prescription of gi- 
ant planet formati on by core accretion presented by 
iMarlcv et al.l (12007H . In their implementation, gas ac¬ 
creting onto the core during the phase of runaway gas ac¬ 
cretion (occurring after the gaseous envelope has reached 
a mass roughly equal to that of the core and has hydro- 
dynamically collapsed onto the core) is able to quickly 
radiate away most of its initial gravitational potential en¬ 
ergy, lowering its specific entropy and causing the gas to 
reach thermal equilibrium with the local radiation field. 
This is the key reason for the significantly lower lumi¬ 
nosity (at early ages) of the “cold-start” models rela¬ 
tive to “hot-start” models. The “hot-start” models of 
iBaraffe et al.l (I2003H require newly-formed giant planets 
to slowly radiate all of their initial gravitational poten¬ 
tial energy, and thus have initially high entropy, as would 
be expec ted for planets formed dir ectly by gravitational 
collapse dSoiegel fc Burrows! 120121 ). For a planet of a 
given luminosity, this amounts to a basic degeneracy in 
the models between its mass, initial conditions, and age, 
such that even if we knew the exact age, it would not 
be clear which model best describes the planet’s evolu¬ 
tionary state since we do not know the exact formation 
mechanism. 

Not only are the luminosities of planets of a given 
mass and age quite disc repant between the models of 
iFortnev et all (1200811 and IBaraffe et al.l (I2003H up to ages 
of some tens of Myr, or even ~ 100 Myr for more mas¬ 
sive planets with m p > 5 Mj up , due to various assump¬ 
tions about their initial conditions, but their spectra 
(and thus their apparent brightness in the various pho¬ 
tometric bands) are also sensitive to the assumed atmo¬ 
spheric conditions (e.g. clouds, composition). Indeed, 
the treatment of atmospheric che mistry and condensa ¬ 
tion differs between the models of lFortnev et all (|2008H 
and lBaraffe et al.l (|2003H , resulting in non-convergence at 
late times when the initial conditions have been forgot¬ 
ten. 

Our assumption that the ages of the planets are the 
same as that of their host could also be a source of un¬ 
certainty in our analysis. Giant planet formation by core 
accretion is believe d to take place over the co urse of sev¬ 
eral Myr (see e.g. iKennedv fc Kenvonl 120081 and refer¬ 
ences therein) before dispersal of the gaseous disk, while 
formation by direct gravitational collapse can occ ur as 
quick ly as a few orbital periods (~ dynamical time: IBossI 
1200011 and up to some ~ 0.5 Myr dHelled et al.l l200(§[ 7 if 


captured planetesimals heat up the contracting proto¬ 
planet. Thus, without knowledge of the exact formation 
mechanism of each planet, we could be overestimating 
its age by up to several Myr when using the “cold-start” 
models. It could also be the case that there are two chan¬ 
nels of planet formation, with some planets being formed 
by core accretion and others by fragmentation. Our hope 
is that by performing our analysis twice, once with each 
class of models, we can bracket the range of expected 
luminosities for a given planet mass and age. 

We also note that these planet evolutionary models 
are poorly calibrated, and comparisons with the few ex¬ 
isting benchmark brown dwarfs have revealed significant 
discrepancies between the mo del predicted luminosities 
and the actual m easurements (|Dupuv et al.l I2009L12014 
iCrem) et al.ll2012t l . These studies have found that evolu¬ 
tionary models systematically underpredict luminosities 
by factors of up to 2 — 3, corresponding to mass overes¬ 
timates between 15 — 30 percent. This caveat should be 
kept in mind when interpreting our conclusions. 

Nevertheless, we do not find a significant difference 
between our results derived from the “hot-” and “cold- 
start” calculations (see S I5.2l for discussion), so we do not 
consider the use (or choice) of an evolutionary model a 
dominant source of uncertainty in our analysis (at least 
not over the region of parameter space we are constrain¬ 
ing). 

6.5. Separation Distribution Function in Microlensing 

Surveys 

The results repo rted by the microlens ing surveys of 
iGould et al.1 (120101) and ISumi et al.1 (l 2010 fl were obtained 
under the assumption that there are equal numbers of 
planets per logarithmic interval in projected separation, 
i.e. dN p i/d logs is a constant. The semimajor axis distri¬ 
bution function we adopt, dN p \/d\oga oc , is related, 
but not equivalent, to d/V p i/<ilogs due to conflation with 
the distributions of lens masses and distances, and pos¬ 
sibly the distribution of projection angles (a function of 
e, i, oj , and Mo). Since ft is one of the parameters in our 
planet population model for whi ch we are deriving con¬ 
straints based on the results of IGould et al.1 (120101) and 
ISumi et aD (|2010f) , their assumption of a log-uniform dis¬ 
tribution in projected separation potentially introduces 
uncertainties for which we do not account. 

To investigate this issue, we model the distributions of 
logs for our planet populations over the range | logs| < 
(l/3)(loggo + 4.3), where q o = 5 x 1 0~ 4 (the appr oxi¬ 
mate range of log s to which the IGould et al.1 L2XIL0I and 
ISumi et all 120101 surveys are sensitive; see § m , as 
d/Vpi/dlogs oc s 7 and measure the value of 7 with a 
least squares fit across the 68 % confidence interval we 
derive on /3 for both the “hot-” and “cold-start” results 
(separately). We plot these best-fit values of 7 as a func¬ 
tion of /3 in figure [301 We find in both the “hot-” and 
“cold-start” cases that our results are c onsistent w ith an 
Opik’s Law distribution (i.e. 7 = 0 : [Polkl fl924l) over 
some range of the 68 % confidence interval we derive on 
/3 (the range of /3 plotted in each panel), although 7 
can differ significantly across the full confidence inte rval. 
Furthermore, the results of the IGould et al.l (I 2010 D and 
ISumi et aD (120101) surveys are only weakly dependent on 
their assumption of Opik’s Law, since their surveys are 
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Figure 30. The slope of the distribution of s measured over the range of s to which the I Gould et al.l (120101 ') and ISumi et all (12010 ) 
microlensing surveys are sensitive, 7 , defined implicitly as dN p \/d log s oc s 7 , as a function of 0. The range of 0 shown in these panels is 
the lcr uncertainty (i.e. the 68 % confidence interval) that we derive from our joint comparison with microlensing, RV, and direct imaging 
surveys. In each panel, the solid, black line shows the values of 7 we find when a, A, and a ou t are held fixed at their median values, and the 
grey band shows the “la range” that we find when a, A, and a ou t are each allowed to vary (independently) within their 68 % confidence 
intervals (see table [2}. The dotted, red, vertical lines represent the median values of 0 we derive from our joint comparison with all surveys, 
and the solid, blue , hor izontal line at 7 = 0 corresponds to a log-uniform distribution of s (i.e. Opik’s Law; Opik 1924), as assumed by 
I Gould et all (2010) and ISumi et all 120114 ). 


sensitive to a relatively small range of s (and thus, a), es¬ 
pecially compared to the range we explore in this study. 
Thus, we argue that while this could be a non-trivial 
source of unquantified uncertainty in our analysis, is it 
probably not a dominant source. 

6 .6. Comparisons with “Intermediate Data Products” 

The most correct and robust way to perform this study 
would be to compare our planet distribution model di¬ 
rectly to the data collected by each respective survey. 
For example, to compare with l ong-term trends dete cted 
by the CPS TRENDS survey of lMontet et alJ (|2014D . we 
would need to follow the procedure outlined below. 

1. For a given planet population, map each ( m P7 a ) 
pair into the observables ( K : P), marginalizing over 
the appropriate stellar and orbital parameters. 

2. Using the derived ( K,P ), generate simulated RV 
curves with the same number of observations and 
cadence for each star in the actual CPS sample, in¬ 
cluding uncertainties due to systematics and stel¬ 
lar jitter. Variations in orbital parameters (includ¬ 
ing phases) are included as part of the mapping 
(m p , a) -A ( K,P). 

3. Pass these simulated RV curves through the detec¬ 
tion pipeline of the CPS TRENDS survey to deter¬ 
mine which planets produce detectable trends. 

4. Compute the probability that the detected trends 
have accelerations consistent with the actual trends 
reported by the CPS. For any trends we expect 
to be deteced by the CPS that are not actually 
detected, we would compute the associated Poisson 
probability. 


5. Multiplying together these probabilities, assign the 
planet population model under consideration a sta¬ 
tistical weight. Repeat this procedure for all pop¬ 
ulations to constrain the parameters of our model, 

{cr, /?, .A, n 0 ut}- 

Similar types of analyses would be performed to com¬ 
pare our planet models to microlensing and direct imag¬ 
ing surveys. Such a study would require all original data 
and detection software for each survey, and thus would 
be much more involved than the methods we have de¬ 
veloped. However, given the levels of uncertainty due 
to the sources discussed previously in this section, we 
believe such a careful analysis is currently unwarranted. 

We instead compare our planet population model to 
“intermediate data products” of the various surveys we 
consider. This results in sources of unquantified uncer¬ 
tainty that are only somewhat alleviated by the fact that 
we use the appropriate detection limits (or estimates 
thereof). Returning to the case of comparison with the 
trends detected by the CPS, our methodology is to use 
estimated detection limits (taking into account the ac¬ 
tual number of observations, stellar jitter, measurement 
uncertainties, and time baseline of observations) for each 
star in their sample in the following manner. The derived 
observables, ( K,P ), for each planet in a given popula¬ 
tion is determined to be either detectable or undetectable 
given the estimated sensitivites. We then assign a statis¬ 
tical weight to the population under consideration based 
on the total number of trends we expect to be detected 
by the CPS TRENDS survey and the likelihood function 
we approximate for their survey (see § [T]) on the total 
(actual) number of detected trends. The “intermediate 
data product” to which we compare for the CPS is thus 
the total number of detected trends. 

7. SUMMARY AND CONCLUSIONS 
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In this paper, we assume a planet distribution func¬ 
tion described by a simple, power-law form with 
four free parameters: d 2 N p \/(d\ogm p dloga) = 
A(m p /Msa.t) a (a/2.5 AU)* 3 , with an outer cutoff radius 
of the separation distribution function of a ou t- We gen¬ 
erate random planet populations from this distribution 
function, map their properties into the observables rel¬ 
evant to microlensing, RV, and direct imaging surveys, 
and compare the expected observations (determined us¬ 
ing the actual detection limits) to the reported results 
of the five surveys we consider. We find that each sur¬ 
vey individually can only constrain combinations of the 
parameters {a, j3, A , a out }, but by performing a simulta¬ 
neous analysis with the results of all surveys, we demon¬ 
strate that there is a single planet population that is 
consistent with them all. 

Each survey is sensitive to different, yet complemen¬ 
tary, regions of parameter space with varying amounts 
of mutual overlap. This is important, as it allows each 
survey to sample the planet distribution function over 
different ranges of planet masses and semimajor axes, 
providing powerful constraints on the properties of plan¬ 
ets at large sepa r ation s. W e find that the com bination of 
the I Gould et all (120101) and lSumi et al.1 (120101) microlens¬ 
ing surveys tightly constrain a and A. With these two 
parameters (essentially) fi xed, i ncluding the results of 
the CPS TRENDS survey dMontet et al.l 120141 ) places a 
lower bound on the allowed values of (3 and a ou t in order 
to exp l ain th e four significant long-term RV drifts that 
iMontet et al.l (120141) (convincingly) argue must be due 
to giant planets. However, we find tha t the subse q uent 
inclu sion of the res ults of the GDPS (lLafreniere et al.1 
120071) and PALMS (jBowler et al.l 120151 ) direct imaging 
surveys can only place relatively weak constraints on the 
values of P and a ou t, and we are left with a large range 
of allowed values on these two parameters. 

The median values and 68% confidence intervals of the 
parameters of a planet population consistent with mi¬ 
crolensing, RV, and direct imaging surveys that we derive 
when assuming “hot-start” planet evolutionary models 
are a = -0.86±gug, P = A = 0.21+g;^ dex“ 2 , 

and Gout = lOj^ y AU, and when assuming “cold-start” 
evolutionary models, we find a = — 0.85to;ig, P = 

l-l±i:|, A = 0.2llo:i5 dex- 2 , and a ou t = 12t^° 2 AU. 

We find that while the results from our analysis as¬ 
suming “hot-start” evolutionary models are slightly more 
constraining than those we find when assuming “cold- 
start” models, the difference is not significant. The lack 
of a significant difference can be attributed to a few rea¬ 
sons: 1) the difference in the predicted luminosities be¬ 
tween the two classes of models decreases with planet 
mass at all ages, 2) our assumed power-law form for the 
planet distribution function, and 3) the most likely range 
of a ou t we find lies within a region where direct imaging 
surveys are mostly insensitive to planets. We find that a 
is relatively well-constrained to negative values, meaning 
the mass function is heavily weighted towards lower-mass 
planets. The predicted luminosities for lower-mass plan¬ 
ets differ less between the “hot-” and “cold-start” mod- 
els relative to the m ore massive planets (see figure 1 of 
iFortnev et al.ll2008 1. Furthermore, because our assumed 
planet distribution function has a power-law form, small 


changes in cc, /?, and/or a ou t can have a dramatic effect on 
the number of expected planet detections. Additionally, 
as the projected separation of a planet approaches the 
inner working angle, its detectability by direct imaging 
surveys dramatically decreases, meaning that we do not 
predict any detections of planets with projected separa¬ 
tions rj_ < 10 AU for either imaging survey, regardless 
of our choice of evolutionary model. 

The degeneracy between /3 and a ou t suffered by our 
comparison with the results of direct imaging surveys, 
where arbitrarily large values of /3 are allowed when 
o 0 ut < 10 AU (where these surveys are insensitive to 
planet detections) and large values of a ou t are allowed 
for small values of /3 (see figure l25l) , limits our ability 
to derive tight constraints on these two parameters. If 
we knew the value of either /3 or a ou tj we would im¬ 
mediately be able to break this degeneracy and derive 
much better constraints. Observational studies mea¬ 
suring (sub-)millimeter-wave continuum emission from 
protoplanetary disks around young M stars have in¬ 
ferred a wide range of characteristic radii, beyond which 
the dust surface density declines exponentia lly , between 
~ 20 — 200 AIJ (1 Andre ws fc William^ l2007t llsell a et al.1 
120091 : [Andrews et al.H2009LI2010t IGuilloteau et al.H201lD . 

The median value and 68% confidence interval we de¬ 
rive on Gout are consistent with these observed disk sizes, 
however the long tail towards larger a ou t where direct 
imaging surveys lose sensitivity due to a limited FOV 
thus seems less likely than we estimate for a ou t > 200 
(see figure 1271) . While planets formed in the disk could be 
placed into such very long-period orbits due to dynamical 
processes, such as planet-planet scatteri ng, they are u n- 

S to end up on stable orbits dDodson-Robinson et al.l 
Futhermore, in a cluster environment, the 
timescale for planet “ionization” by passing stars is 
~ 200 Myr(n/pc _3 ) _1 (a/200 AU) -i (u/0.3 km s -1 ) -1 , 
where n is the number density of stars and v is the rela¬ 
tive velo city between the p l anetary syst em and the pass¬ 
ing star (jAntognini fc Thomr>sonll2015D . 

In this paper, we find that a single, continuous, joint- 
power law distribution function in both planet mass 
and semimajor axis can simultaneously explain existing 
observations from microlensing, RV, and direct imag¬ 
ing. However, we acknowledge the fact that the true 
planet distribution function is probably more compli¬ 
cated, with breaks in both ma ss and semimajor axis, 
and m ultiple power-law slopes. iDodson-Robinson et al.l 
(2009) demonstrate that giant planet formation by gravi¬ 
tational instability can succeed at large separations, but, 
even under the most favorable conditions, in situ giant 
planet formation by core accretion cannot occur beyond 
a ~ 20 AU around M stars. Furthermore, outward mi¬ 
gration of giant planets from the inner disk to larg e sep¬ 
arations seems implausible (iPcpliriski et al.l 120081 ), and 
although planet-planet scattering can shuffle giant plan¬ 
ets out to lar ge separations, it i s unlikely they end up on 
stable orbits dDodson-Robinson et al.ll2009i) . Thus, if the 
mode of giant planet formation switches from core accre¬ 
tion to disk fragmentation as separation from the host 
star increases, it is reasonable to expect that the dis¬ 
tribution of semimajor axes will have at least one break. 
The different formation mechanisms probably lead to dis¬ 
tinct mass functions as well. However, since the existing 
observations can be explained by our simple distribution 
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function, we argue that there is not yet a need to in¬ 
troduce more complicated distribution functions. As di¬ 
rect imaging surveys become more sensitive (e.g. JWST, 
WFIRST-AFTA ), we may find evidence of such breaks 
in the planet distribution function. 

Finally, while we have shown that wide-separation 
(a > 2 AU) planets around M dwarfs can be explained by 
a single population, it would be interesting to see if such 
a population can explain the over-abundan ce of short- 
timesc ale microlensing events identified by iSumi et~ahl 
(2011). These short-timescale events are consistent with 
planetary-mass objects that are either widely separated 
from their hosts, or are freely floating. In a future paper, 
we will investigate the question of whether or not our 
population of wide-separation planets could be responsi¬ 
ble for these short-timescale events, and if not, determine 
what would be required to explain them. 
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